Ticket #8988: trac_8988_add_support_for_toric_varieties.patch

File trac_8988_add_support_for_toric_varieties.patch, 72.6 KB (added by novoselt, 12 years ago)

Apply this patch and doctest fix.

  • doc/en/reference/schemes.rst

    # HG changeset patch
    # User Andrey Novoseltsev <novoselt@gmail.com>
    # Date 1274174310 21600
    # Node ID 5c1acb93c35758e535485d50448299f3623e7c84
    # Parent  e35a4895803ac221797c1a51494234cd4675cb86
    Trac 8988: Add support for toric varieties.
    
    diff -r e35a4895803a -r 5c1acb93c357 doc/en/reference/schemes.rst
    a b  
    1616   sage/schemes/generic/ambient_space
    1717   sage/schemes/generic/affine_space
    1818   sage/schemes/generic/projective_space
     19   sage/schemes/generic/toric_variety
    1920   sage/schemes/generic/algebraic_scheme
    2021   sage/schemes/generic/hypersurface
    2122
  • sage/geometry/fan.py

    diff -r e35a4895803a -r 5c1acb93c357 sage/geometry/fan.py
    a b  
    770770
    771771        .. NOTE::
    772772
    773             "Direct call" syntax is a synonym for :meth:`cones` method.
     773            "Direct call" syntax is a synonym for :meth:`cones` method except
     774            that in the case of no input parameters this function returns
     775            just ``self``.
    774776
    775777        INPUT:
    776778
     
    780782
    781783        OUTPUT:
    782784
    783         - cones of ``self`` of the specified (co)dimension.
     785        - cones of ``self`` of the specified (co)dimension if it was given,
     786          otherwise ``self``.
    784787
    785788        TESTS::
    786789
     
    802805            ...
    803806            ValueError: dimension and codimension
    804807            cannot be specified together!
     808            sage: fan() is fan
     809            True
    805810        """
    806         return self.cones(dim, codim)
     811        if dim is None and codim is None:
     812            # "self.cones()" returns all cones, but for the call syntax
     813            # "self()" we return just "self", which seems to be more natural
     814            # and convenient for ToricVariety.fan() method.
     815            return self
     816        else:
     817            return self.cones(dim, codim)
    807818
    808819    def __cmp__(self, right):
    809820        r"""
  • sage/schemes/generic/all.py

    diff -r e35a4895803a -r 5c1acb93c357 sage/schemes/generic/all.py
    a b  
    77from morphism         import is_SchemeMorphism
    88from projective_space import ProjectiveSpace, is_ProjectiveSpace
    99from scheme           import is_Scheme, is_AffineScheme
     10from toric_variety import AffineToricVariety, ToricVariety
    1011from hypersurface     import ProjectiveHypersurface, AffineHypersurface
    1112
    1213
  • new file sage/schemes/generic/toric_variety.py

    diff -r e35a4895803a -r 5c1acb93c357 sage/schemes/generic/toric_variety.py
    - +  
     1r"""
     2Toric varieties
     3
     4This module provides support for (normal) toric varieties, corresponding to
     5:class:`rational polyhedral fans <sage.geometry.fan.RationalPolyhedralFan>`.
     6See also :mod:`~sage.schemes.generic.fano_toric_variety` for a more
     7restrictive class of (weak) Fano toric varieties.
     8
     9An **excellent reference on toric varieties** is the book "Toric Varieties" by
     10David A. Cox, John B. Little, and Hal Schenck [CLS]_. Its draft **is freely
     11available** at
     12http://www.cs.amherst.edu/~dac/toric.html
     13**but will be removed** from this site once it is published, so hurry up!
     14
     15The interface to this module is provided through functions
     16:func:`AffineToricVariety` and :func:`ToricVariety`, although you may also be
     17interested in :func:`normalize_names`.
     18
     19.. NOTE::
     20
     21    We do NOT build "general toric varieties" from affine toric varieties.
     22    Instead, we are using the quotient representation of toric varieties with
     23    the homogeneous coordinate ring (a.k.a. Cox's ring or the total coordinate
     24    ring). This description works best for simplicial fans of the full
     25    dimension.
     26
     27REFERENCES:
     28
     29..  [CLS]
     30    David A. Cox, John B. Little,  Hal Schenck,
     31    "Toric Varieties", Graduate Studies in Mathematics,
     32    Amer. Math. Soc., Providence, RI, to appear.
     33
     34AUTHORS:
     35
     36- Andrey Novoseltsev (2010-05-17): initial version.
     37
     38EXAMPLES:
     39
     40We start with constructing the affine plane as an affine toric variety. First,
     41we need to have a corresponding cone::
     42
     43    sage: quadrant = Cone([(1,0), (0,1)])
     44
     45If you don't care about variable names and the base field, that's all we need
     46for now::
     47
     48    sage: A2 = AffineToricVariety(quadrant)
     49    sage: A2
     50    2-d affine toric variety
     51    sage: origin = A2(0,0)
     52    sage: origin
     53    [0 : 0]
     54
     55Only affine toric varieties have points whose (homogeneous) coordinates
     56are all zero. ::
     57
     58    sage: parent(origin)
     59    Set of Rational Points of 2-d affine toric variety
     60
     61As you can see, by default toric varieties live over the field of rational
     62numbers::
     63
     64    sage: A2.base_ring()
     65    Rational Field
     66
     67While usually toric varieties are considered over the field of complex
     68numbers, for computational purposes it is more convenient to work with fields
     69that have exact representation on computers. You can also always do ::
     70
     71    sage: C2 = AffineToricVariety(quadrant, base_field=CC)
     72    sage: C2.base_ring()
     73    Complex Field with 53 bits of precision
     74    sage: C2(1,2+i)
     75    [1.00000000000000 : 2.00000000000000 + 1.00000000000000*I]
     76
     77or even ::
     78
     79    sage: F = CC["a, b"].fraction_field()
     80    sage: F.inject_variables()
     81    Defining a, b
     82    sage: A2 = AffineToricVariety(quadrant, base_field=F)
     83    sage: A2(a,b)
     84    [a : b]
     85
     86OK, if you need to work only with affine spaces,
     87:func:`~sage.schemes.generic.affine_space.AffineSpace` may be a better way to
     88construct them. Our next example is the product of two projective lines
     89realized as the toric variety associated to the
     90:func:`face fan <sage.geometry.fan.FaceFan>` of the "diamond"::
     91
     92    sage: diamond = lattice_polytope.octahedron(2)
     93    sage: diamond.vertices()
     94    [ 1  0 -1  0]
     95    [ 0  1  0 -1]
     96    sage: fan = FaceFan(diamond)
     97    sage: P1xP1 = ToricVariety(fan)
     98    sage: P1xP1
     99    2-d toric variety covered by 4 affine patches
     100    sage: P1xP1.fan().ray_matrix()
     101    [ 1  0 -1  0]
     102    [ 0  1  0 -1]
     103    sage: P1xP1.gens()
     104    (z0, z1, z2, z3)
     105
     106We got four coordinates - two for each of the projective lines, but their
     107names are perhaps not very well chosen. Let's make `(x,y)` to be coordinates
     108on the first line and `(s,t)` on the second one::
     109
     110    sage: P1xP1 = ToricVariety(fan, coordinate_names="x s y t")
     111    sage: P1xP1.gens()
     112    (x, s, y, t)
     113
     114Now, if we want to define subschemes of this variety, the defining polynomials
     115must be homogeneous in each of these pairs::
     116
     117    sage: P1xP1.inject_variables()
     118    Defining x, s, y, t
     119    sage: P1xP1.subscheme(x)
     120    Closed subscheme of 2-d toric variety
     121    covered by 4 affine patches defined by:
     122      x
     123    sage: P1xP1.subscheme(x^2 + y^2)
     124    Closed subscheme of 2-d toric variety
     125    covered by 4 affine patches defined by:
     126      x^2 + y^2
     127    sage: P1xP1.subscheme(x^2 + s^2)
     128    Traceback (most recent call last):
     129    ...
     130    ValueError: x^2 + s^2 is not homogeneous
     131    on 2-d toric variety covered by 4 affine patches!
     132    sage: P1xP1.subscheme([x^2*s^2 + x*y*t^2 +y^2*t^2, s^3 + t^3])
     133    Closed subscheme of 2-d toric variety
     134    covered by 4 affine patches defined by:
     135      x^2*s^2 + x*y*t^2 + y^2*t^2,
     136      s^3 + t^3
     137
     138While we don't build toric varieties from affine toric varieties, we still can
     139access the "building pieces"::
     140
     141    sage: patch = P1xP1.affine_patch(2)
     142    sage: patch
     143    2-d affine toric variety
     144    sage: patch.fan().ray_matrix()
     145    [1 0]
     146    [0 1]
     147    sage: patch.embedding_morphism()
     148    Scheme morphism:
     149      From: 2-d affine toric variety
     150      To:   2-d toric variety covered by 4 affine patches
     151      Defn: Defined on coordinates by sending [x : s] to
     152            [x : s : 1 : 1]
     153
     154The patch above was specifically chosen to coincide with our representation of
     155the affine plane before, but you can get the other three patches as well.
     156(While any cone of a fan will correspond to an affine toric variety, the main
     157interest is usually in the generating fans as "the biggest" affine
     158subvarieties, and these are precisely the patches that you can get from
     159:meth:`~ToricVariety_field.affine_patch`.)
     160
     161All two-dimensional toric varieties are "quite nice" because any
     162two-dimensional cone is generated by exactly two rays. From the point of view
     163of the corresponding toric varieties, this means that they have at worst
     164quotient singularities::
     165
     166    sage: P1xP1.is_orbifold()
     167    True
     168    sage: P1xP1.is_smooth()
     169    True
     170    sage: TV = ToricVariety(NormalFan(diamond))
     171    sage: TV.fan().ray_matrix()
     172    [-1  1 -1  1]
     173    [ 1  1 -1 -1]
     174    sage: TV.is_orbifold()
     175    True
     176    sage: TV.is_smooth()
     177    False
     178
     179In higher dimensions worse things can happen::
     180
     181    sage: TV3 = ToricVariety(NormalFan(lattice_polytope.octahedron(3)))
     182    sage: TV3.fan().ray_matrix()
     183    [-1  1 -1  1 -1  1 -1  1]
     184    [-1 -1  1  1 -1 -1  1  1]
     185    [ 1  1  1  1 -1 -1 -1 -1]
     186    sage: TV3.is_orbifold()
     187    False
     188
     189Fortunately, we can perform a (partial) resolution::
     190
     191    sage: TV3_res = TV3.resolve_to_orbifold()
     192    sage: TV3_res.is_orbifold()
     193    True
     194    sage: TV3_res.fan().ngenerating_cones()
     195    12
     196    sage: TV3.fan().ngenerating_cones()
     197    6
     198
     199In this example we had to double the number of affine patches. The result is
     200still singular::
     201
     202    sage: TV3_res.is_smooth()
     203    False
     204
     205You can resolve it further using :meth:`~ToricVariety_field.resolve` method,
     206but (at least for now) you will have to specify which rays should be inserted
     207into the fan. See also
     208:func:`~sage.schemes.generic.fano_toric_variety.FanoToricVariety`,
     209which can construct some other "nice partial resolutions."
     210
     211Below you will find detailed descriptions of available functions. If you are
     212familiar with toric geometry, you will likely see that many important objects
     213and operations are unavailable. However, this module is under active
     214development and hopefully will improve in future releases of Sage. If there
     215are some particular features that you would like to see implemented ASAP,
     216please consider reporting them to the Sage Development Team or even
     217implementing them on your own as a patch for inclusion!
     218"""
     219
     220
     221#*****************************************************************************
     222#       Copyright (C) 2010 Andrey Novoseltsev <novoselt@gmail.com>
     223#       Copyright (C) 2010 William Stein <wstein@gmail.com>
     224#
     225#  Distributed under the terms of the GNU General Public License (GPL)
     226#
     227#                  http://www.gnu.org/licenses/
     228#*****************************************************************************
     229
     230
     231from sage.geometry.all import Cone, Fan
     232from sage.matrix.all import matrix
     233from sage.misc.all import latex
     234from sage.modules.free_module_element import vector
     235from sage.rings.all import PolynomialRing, QQ, is_Field, is_FractionField
     236from sage.schemes.generic.algebraic_scheme import AlgebraicScheme_subscheme
     237from sage.schemes.generic.ambient_space import AmbientSpace
     238from sage.schemes.generic.homset import SchemeHomset_coordinates
     239from sage.schemes.generic.morphism import (SchemeMorphism_coordinates,
     240                                           SchemeMorphism_on_points,
     241                                           is_SchemeMorphism)
     242from sage.schemes.generic.scheme import is_Scheme
     243from sage.structure.sequence import Sequence
     244
     245
     246# Default prefix for indexed coordinates
     247DEFAULT_PREFIX = "z"
     248
     249
     250def is_ToricVariety(x):
     251    r"""
     252    Check if ``x`` is a toric variety.
     253
     254    INPUT:
     255
     256    - ``x`` -- anything.
     257
     258    OUTPUT:
     259
     260    - ``True`` if ``x`` is a :class:`toric variety <ToricVariety_field>` and
     261      ``False`` otherwise.
     262
     263    .. NOTE::
     264
     265        While projective spaces are toric varieties mathematically, they are
     266        not toric varieties in Sage due to efficiency considerations, so this
     267        function will return ``False``.
     268
     269    EXAMPLES::
     270
     271        sage: from sage.schemes.generic.toric_variety import is_ToricVariety
     272        sage: is_ToricVariety(1)
     273        False
     274        sage: fan = FaceFan(lattice_polytope.octahedron(2))
     275        sage: P = ToricVariety(fan)
     276        sage: P
     277        2-d toric variety covered by 4 affine patches
     278        sage: is_ToricVariety(P)
     279        True
     280        sage: is_ToricVariety(ProjectiveSpace(2))
     281        False
     282    """
     283    return isinstance(x, ToricVariety_field)
     284
     285
     286def ToricVariety(fan,
     287                 coordinate_names=None,
     288                 coordinate_indices=None,
     289                 base_field=QQ):
     290    r"""
     291    Construct a toric variety.
     292
     293    INPUT:
     294
     295    - ``fan`` -- :class:`rational polyhedral fan
     296      <sage.geometry.fan.RationalPolyhedralFan>`;
     297
     298    - ``coordinate_names`` -- names of variables for the coordinate ring, see
     299      :func:`normalize_names` for acceptable formats. If not given, indexed
     300      variable names will be created automatically;
     301
     302    - ``coordinate_indices`` -- list of integers, indices for indexed
     303      variables. If not given, the index of each variable will coincide with
     304      the index of the corresponding ray of the fan;
     305
     306    - ``base_field`` -- base field of the toric variety (default: `\QQ`).
     307
     308    OUTPUT:
     309
     310    - :class:`toric variety <ToricVariety_field>`.
     311
     312    EXAMPLES:
     313
     314    We will create the product of two projective lines::
     315
     316        sage: fan = FaceFan(lattice_polytope.octahedron(2))
     317        sage: fan.ray_matrix()
     318        [ 1  0 -1  0]
     319        [ 0  1  0 -1]
     320        sage: P1xP1 = ToricVariety(fan)
     321        sage: P1xP1.gens()
     322        (z0, z1, z2, z3)
     323
     324    Let's create some points::
     325
     326        sage: P1xP1(1,1,1,1)
     327        [1 : 1 : 1 : 1]
     328        sage: P1xP1(0,1,1,1)
     329        [0 : 1 : 1 : 1]
     330        sage: P1xP1(0,1,0,1)
     331        Traceback (most recent call last):
     332        ...
     333        TypeError: coordinates (0, 1, 0, 1)
     334        are in the exceptional set!
     335
     336    We cannot set to zero both coordinates of the same projective line!
     337
     338    Let's change the names of the variables. We have to re-create our toric
     339    variety::
     340
     341        sage: P1xP1 = ToricVariety(fan, "x s y t")
     342        sage: P1xP1.gens()
     343        (x, s, y, t)
     344
     345    Now `(x, y)` correspond to one line and `(s, t)` to the other one. ::
     346
     347        sage: P1xP1.inject_variables()
     348        Defining x, s, y, t
     349        sage: P1xP1.subscheme(x*s-y*t)
     350        Closed subscheme of 2-d toric variety
     351        covered by 4 affine patches defined by:
     352          x*s - y*t
     353    """
     354    if not is_Field(base_field):
     355        raise TypeError("need a field to construct a toric variety!\n Got %s"
     356                        % base_field)
     357    return ToricVariety_field(fan, coordinate_names, coordinate_indices,
     358                             base_field)
     359
     360
     361def AffineToricVariety(cone, *args, **kwds):
     362    r"""
     363    Construct an affine toric variety.
     364
     365    INPUT:
     366
     367    - ``cone`` -- :class:`strictly convex rational polyhedral cone
     368      <sage.geometry.cone.ConvexRationalPolyhedralCone>`.
     369
     370    This cone will be used to construct a :class:`rational polyhedral fan
     371    <sage.geometry.fan.RationalPolyhedralFan>`, which will be passed to
     372    :func:`ToricVariety` with the rest of positional and keyword arguments.
     373
     374    OUTPUT:
     375
     376    - :class:`toric variety <ToricVariety_field>`.
     377
     378    .. NOTE::
     379
     380        The generating rays of the fan of this variety are guaranteed to be
     381        listed in the same order as the rays of the original cone.
     382
     383    EXAMPLES:
     384
     385    We will create the affine plane as an affine toric variety::
     386
     387        sage: quadrant = Cone([(1,0), (0,1)])
     388        sage: A2 = AffineToricVariety(quadrant)
     389        sage: origin = A2(0,0)
     390        sage: origin
     391        [0 : 0]
     392        sage: parent(origin)
     393        Set of Rational Points of 2-d affine toric variety
     394
     395    Only affine toric varieties have points whose (homogeneous) coordinates
     396    are all zero.
     397    """
     398    if not cone.is_strictly_convex():
     399        raise ValueError("affine toric varieties are defined for strictly "
     400                         "convex cones only!")
     401    # We make sure that Fan constructor does not meddle with the order of
     402    # rays, this is very important for affine patches construction
     403    fan = Fan([tuple(range(cone.nrays()))], cone.rays(),
     404              check=False, normalize=False)
     405    return ToricVariety(fan, *args, **kwds)
     406
     407
     408class ToricVariety_field(AmbientSpace):
     409    r"""
     410    Construct a toric variety associated to a rational polyhedral fan.
     411
     412    .. WARNING::
     413
     414        This class does not perform any checks of correctness of input. Use
     415        :func:`ToricVariety` and :func:`AffineToricVariety` to construct toric
     416        varieties.
     417
     418    INPUT:
     419
     420    - ``fan`` -- :class:`rational polyhedral fan
     421      <sage.geometry.fan.RationalPolyhedralFan>`;
     422
     423    - ``coordinate_names`` -- names of variables, see :func:`normalize_names`
     424      for acceptable formats. If ``None``, indexed variable names will be
     425      created automatically;
     426
     427    - ``coordinate_indices`` -- list of integers, indices for indexed
     428      variables. If ``None``, the index of each variable will coincide with
     429      the index of the corresponding ray of the fan;
     430
     431    - ``base_field`` -- base field of the toric variety.
     432
     433    OUTPUT:
     434
     435    - :class:`toric variety <ToricVariety_field>`.
     436
     437    TESTS::
     438
     439        sage: fan = FaceFan(lattice_polytope.octahedron(2))
     440        sage: P1xP1 = ToricVariety(fan)
     441    """
     442
     443    def __init__(self, fan, coordinate_names, coordinate_indices, base_field):
     444        r"""
     445        See :class:`ToricVariety_field` for documentation.
     446
     447        TESTS::
     448
     449            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     450            sage: P1xP1 = ToricVariety(fan)
     451        """
     452        self._fan = fan
     453        super(ToricVariety_field, self).__init__(fan.lattice_dim(),
     454                                                 base_field)
     455        self._torus_factor_dim = fan.lattice_dim() - fan.dim()
     456        coordinate_names = normalize_names(coordinate_names,
     457                        fan.nrays() + self._torus_factor_dim, DEFAULT_PREFIX,
     458                        coordinate_indices, return_prefix=True)
     459        # Save the prefix for use in resolutions
     460        self._coordinate_prefix = coordinate_names.pop()
     461        self._assign_names(names=coordinate_names, normalize=False)
     462
     463    def __cmp__(self, right):
     464        r"""
     465        Compare ``self`` and ``right``.
     466
     467        INPUT:
     468
     469        - ``right`` -- anything.
     470
     471        OUTPUT:
     472
     473        - 0 if ``right`` is of the same type as ``self``, their fans are the
     474          same, names of variables are the same and stored in the same order,
     475          and base fields are the same. 1 or -1 otherwise.
     476
     477        TESTS::
     478
     479            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     480            sage: P1xP1 = ToricVariety(fan)
     481            sage: P1xP1a = ToricVariety(fan, "x s y t")
     482            sage: P1xP1b = ToricVariety(fan)
     483            sage: cmp(P1xP1, P1xP1a)
     484            1
     485            sage: cmp(P1xP1a, P1xP1)
     486            -1
     487            sage: cmp(P1xP1, P1xP1b)
     488            0
     489            sage: P1xP1 is P1xP1b
     490            False
     491            sage: cmp(P1xP1, 1)
     492            -1
     493        """
     494        c = cmp(type(self), type(right))
     495        if c:
     496            return c
     497        return cmp([self.fan(),
     498                    self.variable_names(),
     499                    self.base_ring()],
     500                   [right.fan(),
     501                    right.variable_names(),
     502                    right.base_ring()])
     503
     504    def _an_element_(self):
     505        r"""
     506        Construct an element of ``self``.
     507
     508        This function is needed (in particular) for the test framework.
     509
     510        OUTPUT:
     511
     512        - a point of ``self`` with coordinates [1 : 2: ... : n].
     513
     514        TESTS::
     515
     516            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     517            sage: P1xP1 = ToricVariety(fan)
     518            sage: P1xP1._an_element_()
     519            [1 : 2 : 3 : 4]
     520        """
     521        return self(range(1, self.ngens() + 1))
     522
     523    def _check_satisfies_equations(self, coordinates):
     524        r"""
     525        Check if ``coordinates`` define a valid point of ``self``.
     526
     527        INPUT:
     528
     529        - ``coordinates`` -- list of elements of the base field of ``self``.
     530
     531        OUTPUT:
     532
     533        - ``True`` if ``coordinates`` do define a valid point of ``self``,
     534          otherwise a ``TypeError`` or ``ValueError`` exception is raised.
     535
     536        TESTS::
     537
     538            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     539            sage: P1xP1 = ToricVariety(fan)
     540            sage: P1xP1._check_satisfies_equations([1,1,1,1])
     541            True
     542            sage: P1xP1._check_satisfies_equations([0,0,1,1])
     543            True
     544            sage: P1xP1._check_satisfies_equations([0,1,0,1])
     545            Traceback (most recent call last):
     546            ...
     547            TypeError: coordinates (0, 1, 0, 1)
     548            are in the exceptional set!
     549            sage: P1xP1._check_satisfies_equations([1,1,1])
     550            Traceback (most recent call last):
     551            ...
     552            TypeError: coordinates (1, 1, 1) must have 4 components!
     553            sage: P1xP1._check_satisfies_equations(1)
     554            Traceback (most recent call last):
     555            ...
     556            TypeError: 1 can not be used as coordinates!
     557            Use a list or a tuple.
     558            sage: P1xP1._check_satisfies_equations([1,1,1,fan])
     559            Traceback (most recent call last):
     560            ...
     561            TypeError: coordinate Rational polyhedral fan
     562            in 2-d lattice N is not an element of Rational Field!
     563        """
     564        try:
     565            coordinates = tuple(coordinates)
     566        except TypeError:
     567            raise TypeError("%s can not be used as coordinates! "
     568                            "Use a list or a tuple." % coordinates)
     569        n = self.ngens()
     570        if len(coordinates) != n:
     571            raise TypeError("coordinates %s must have %d components!"
     572                            % (coordinates, n))
     573        base_field = self.base_ring()
     574        for coordinate in coordinates:
     575            if coordinate not in base_field:
     576                raise TypeError("coordinate %s is not an element of %s!"
     577                                % (coordinate, base_field))
     578        zero_positions = set(position
     579                            for position, coordinate in enumerate(coordinates)
     580                            if coordinate == 0)
     581        if not zero_positions:
     582            return True
     583        for i in range(n - self._torus_factor_dim, n):
     584            if i in zero_positions:
     585                raise ValueError("coordinates on the torus factor cannot be "
     586                                 "zero! Got %s" % str(coordinates))
     587        if len(zero_positions) == 1:
     588            return True
     589        fan = self.fan()
     590        possible_charts = set(fan._ray_to_cones(zero_positions.pop()))
     591        for i in zero_positions:
     592            possible_charts.intersection_update(fan._ray_to_cones(i))
     593        if possible_charts:
     594            return True     # All zeros are inside one generating cone
     595        raise TypeError("coordinates %s are in the exceptional set!"
     596                        % str(coordinates)) # Need str, coordinates is a tuple
     597
     598    def _homset_class(self, *args, **kwds):
     599        r"""
     600        Construct a `Hom`-space for ``self``.
     601
     602        INPUT:
     603
     604        - same as for :class:`SchemeHomset_toric_coordinates_field`.
     605
     606        OUPUT:
     607
     608        - :class:`SchemeHomset_toric_coordinates_field`.
     609
     610        TESTS::
     611
     612            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     613            sage: P1xP1 = ToricVariety(fan)
     614            sage: P1xP1._homset_class(P1xP1, P1xP1.base_ring())
     615            Set of Rational Points of 2-d toric variety
     616            covered by 4 affine patches
     617        """
     618        return SchemeHomset_toric_coordinates_field(*args, **kwds)
     619
     620    def _latex_(self):
     621        r"""
     622        Return a LaTeX representation of ``self``.
     623
     624        OUTPUT:
     625
     626        - string.
     627
     628        TESTS::
     629
     630            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     631            sage: P1xP1 = ToricVariety(fan)
     632            sage: P1xP1._latex_()
     633            '\\mathbb{X}_{\\Sigma^{2}}'
     634        """
     635        return r"\mathbb{X}_{%s}" % latex(self.fan())
     636
     637    def _latex_generic_point(self, coordinates=None):
     638        """
     639        Return a LaTeX representation of a point of ``self``.
     640
     641        INPUT:
     642
     643        - ``coordinates`` -- list of coordinates of a point of ``self``.
     644          If not given, names of coordinates of ``self`` will be used.
     645
     646        OUTPUT:
     647
     648        - string.
     649
     650        EXAMPLES::
     651
     652            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     653            sage: P1xP1 = ToricVariety(fan)
     654            sage: P1xP1._latex_generic_point()
     655            '\\left[z_{0} : z_{1} : z_{2} : z_{3}\\right]'
     656            sage: P1xP1._latex_generic_point([1,2,3,4])
     657            '\\left[1 : 2 : 3 : 4\\right]'
     658        """
     659        if coordinates is None:
     660            coordinates = self.gens()
     661        return r"\left[%s\right]" % (" : ".join(str(latex(coord))
     662                                                for coord in coordinates))
     663
     664    def _point_class(self, *args, **kwds):
     665        r"""
     666        Construct a point of ``self``.
     667
     668        INPUT:
     669
     670        - same as for :class:`SchemeMorphism_toric_coordinates_field`.
     671
     672        OUPUT:
     673
     674        - :class:`SchemeMorphism_toric_coordinates_field`.
     675
     676        TESTS::
     677
     678            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     679            sage: P1xP1 = ToricVariety(fan)
     680            sage: P1xP1._point_class(P1xP1, [1,2,3,4])
     681            [1 : 2 : 3 : 4]
     682        """
     683        return SchemeMorphism_toric_coordinates_field(*args, **kwds)
     684
     685    def _point_morphism_class(self, *args, **kwds):
     686        r"""
     687        Construct a morphism determined by action on points of ``self``.
     688
     689        INPUT:
     690
     691        - same as for :class:`SchemeMorphism_on_points_toric_variety`.
     692
     693        OUPUT:
     694
     695        - :class:`SchemeMorphism_on_points_toric_variety`.
     696
     697        TESTS::
     698
     699            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     700            sage: P1xP1 = ToricVariety(fan)
     701            sage: P1xP1.inject_variables()
     702            Defining z0, z1, z2, z3
     703            sage: P1 = P1xP1.subscheme(z0-z2)
     704            sage: H = P1xP1.Hom(P1)
     705            sage: P1xP1._point_morphism_class(H, [z0,z1,z0,z3])
     706            Scheme morphism:
     707              From: 2-d toric variety covered by 4 affine patches
     708              To:   Closed subscheme of 2-d toric variety
     709              covered by 4 affine patches defined by:
     710              z0 - z2
     711              Defn: Defined on coordinates by sending
     712                    [z0 : z1 : z2 : z3] to [z0 : z1 : z0 : z3]
     713        """
     714        return SchemeMorphism_on_points_toric_variety(*args, **kwds)
     715
     716    def _repr_(self):
     717        r"""
     718        Return a string representation of ``self``.
     719
     720        OUTPUT:
     721
     722        - string.
     723
     724        TESTS::
     725
     726            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     727            sage: P1xP1 = ToricVariety(fan)
     728            sage: P1xP1._repr_()
     729            '2-d toric variety covered by 4 affine patches'
     730        """
     731        result = "%d-d" % self.dimension_relative()
     732        if self.fan().ngenerating_cones() == 1:
     733            result += " affine toric variety"
     734        else:
     735            result += (" toric variety covered by %d affine patches"
     736                       % self.fan().ngenerating_cones())
     737        return result
     738
     739    def _repr_generic_point(self, coordinates=None):
     740        r"""
     741        Return a string representation of a point of ``self``.
     742
     743        INPUT:
     744
     745        - ``coordinates`` -- list of coordinates of a point of ``self``.
     746          If not given, names of coordinates of ``self`` will be used.
     747
     748        OUTPUT:
     749
     750        - string.
     751
     752        EXAMPLES::
     753
     754            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     755            sage: P1xP1 = ToricVariety(fan)
     756            sage: P1xP1._repr_generic_point()
     757            '[z0 : z1 : z2 : z3]'
     758            sage: P1xP1._repr_generic_point([1,2,3,4])
     759            '[1 : 2 : 3 : 4]'
     760        """
     761        if coordinates is None:
     762            coordinates = self.gens()
     763        return "[%s]" % (" : ".join(str(coord) for coord in coordinates))
     764
     765    def _validate(self, polynomials):
     766        """
     767        Check if ``polynomials`` define valid functions on ``self``.
     768
     769        Since this is a toric variety, polynomials must be homogeneous in the
     770        total coordinate ring of ``self``.
     771
     772        INPUT:
     773
     774        - ``polynomials`` -- list of polynomials in the coordinate ring of
     775          ``self`` (this function does not perform any conversions).
     776
     777        OUTPUT:
     778
     779        - ``polynomials`` (the input parameter without any modifications) if
     780          ``polynomials`` do define valid polynomial functions on ``self``,
     781          otherwise a ``ValueError`` exception is raised.
     782
     783        TESTS:
     784
     785        We will use the product of two projective lines with coordinates
     786        `(x, y)` for one and `(s, t)` for the other::
     787
     788            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     789            sage: fan.ray_matrix()
     790            [ 1  0 -1  0]
     791            [ 0  1  0 -1]
     792            sage: P1xP1 = ToricVariety(fan, "x s y t")
     793            sage: P1xP1.inject_variables()
     794            Defining x, s, y, t
     795            sage: P1xP1._validate([x - y, x*s + y*t])
     796            [x - y, x*s + y*t]
     797            sage: P1xP1._validate([x + s])
     798            Traceback (most recent call last):
     799            ...
     800            ValueError: x + s is not homogeneous on
     801            2-d toric variety covered by 4 affine patches!
     802        """
     803        for p in polynomials:
     804            if not self.is_homogeneous(p):
     805                raise ValueError("%s is not homogeneous on %s!" % (p, self))
     806        return polynomials
     807
     808    def affine_patch(self, i):
     809        r"""
     810        Return the ``i``-th affine patch of ``self``.
     811
     812        INPUT:
     813
     814        - ``i`` -- integer, index of a generating cone of the fan of ``self``.
     815
     816        OUTPUT:
     817
     818        - affine :class:`toric variety <ToricVariety_field>` corresponding to
     819          the ``i``-th generating cone of the fan of ``self``.
     820
     821        The result is cached, so the ``i``-th patch is always the same object
     822        in memory.
     823
     824        EXAMPLES::
     825
     826            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     827            sage: P1xP1 = ToricVariety(fan, "x s y t")
     828            sage: patch0 = P1xP1.affine_patch(0)
     829            sage: patch0
     830            2-d affine toric variety
     831            sage: patch0.embedding_morphism()
     832            Scheme morphism:
     833              From: 2-d affine toric variety
     834              To:   2-d toric variety covered by 4 affine patches
     835              Defn: Defined on coordinates by sending [x : t] to
     836                    [x : 1 : 1 : t]
     837            sage: patch1 = P1xP1.affine_patch(1)
     838            sage: patch1.embedding_morphism()
     839            Scheme morphism:
     840              From: 2-d affine toric variety
     841              To:   2-d toric variety covered by 4 affine patches
     842              Defn: Defined on coordinates by sending [y : t] to
     843                    [1 : 1 : y : t]
     844            sage: patch1 is P1xP1.affine_patch(1)
     845            True
     846        """
     847        i = int(i)   # implicit type checking
     848        try:
     849            return self._affine_patches[i]
     850        except AttributeError:
     851            self._affine_patches = dict()
     852        except KeyError:
     853            pass
     854        cone = self.fan().generating_cone(i)
     855        names = self.variable_names()
     856        # Number of "honest fan coordinates"
     857        n = self.fan().nrays()
     858        # Number of "torus factor coordinates"
     859        t = self._torus_factor_dim
     860        names = ([names[ray] for ray in cone.ambient_ray_indices()]
     861                 + list(names[n:]))
     862        patch = AffineToricVariety(cone, names, base_field=self.base_ring())
     863        embedding_coordinates = [1] * n
     864        for k, ray in enumerate(cone.ambient_ray_indices()):
     865            embedding_coordinates[ray] = patch.gen(k)
     866        if t > 0: # Passing "-0" gives unintended result
     867            embedding_coordinates.extend(patch.gens()[-t:])
     868        patch._embedding_morphism = patch.hom(embedding_coordinates, self)
     869        self._affine_patches[i] = patch
     870        return patch
     871
     872    def change_ring(self, F):
     873        r"""
     874        Return a toric variety over ``F`` and otherwise the same as ``self``.
     875
     876        INPUT:
     877
     878        - ``F`` -- field.
     879
     880        OUTPUT:
     881
     882        - :class:`toric variety <ToricVariety_field>` over ``F``.
     883
     884        .. NOTE::
     885
     886            There is no need to have any relation between ``F`` and the base
     887            field of ``self``. If you do want to have such a relation, use
     888            :meth:`base_extend` instead.
     889
     890        EXAMPLES::
     891
     892            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     893            sage: P1xP1 = ToricVariety(fan)
     894            sage: P1xP1.base_ring()
     895            Rational Field
     896            sage: P1xP1_RR = P1xP1.change_ring(RR)
     897            sage: P1xP1_RR.base_ring()
     898            Real Field with 53 bits of precision
     899            sage: P1xP1_QQ = P1xP1_RR.change_ring(QQ)
     900            sage: P1xP1_QQ.base_ring()
     901            Rational Field
     902            sage: P1xP1_RR.base_extend(QQ)
     903            Traceback (most recent call last):
     904            ...
     905            ValueError: no natural map from the base ring
     906            (=Real Field with 53 bits of precision)
     907            to R (=Rational Field)!
     908        """
     909        if self.base_ring() == F:
     910            return self
     911        else:
     912            return ToricVariety(self.fan(), self.variable_names(),
     913                                base_field=F)
     914
     915    def coordinate_ring(self):
     916        r"""
     917        Return the coordinate ring of ``self``.
     918
     919        For toric varieties this is the homogeneous coordinate ring (a.k.a.
     920        Cox's ring and total ring).
     921
     922        OUTPUT:
     923
     924        - polynomial ring.
     925
     926        EXAMPLES::
     927
     928            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     929            sage: P1xP1 = ToricVariety(fan)
     930            sage: P1xP1.coordinate_ring()
     931            Multivariate Polynomial Ring in z0, z1, z2, z3
     932            over Rational Field
     933        """
     934        if "_coordinate_ring" not in self.__dict__:
     935            self._coordinate_ring = PolynomialRing(self.base_ring(),
     936                                                   self.variable_names())
     937        return self._coordinate_ring
     938
     939    def embedding_morphism(self):
     940        r"""
     941        Return the default embedding morphism of ``self``.
     942
     943        Such a morphism is always defined for an affine patch of a toric
     944        variety (which is also a toric varieties itself).
     945
     946        OUTPUT:
     947
     948        - :class:`scheme morphism <SchemeMorphism_on_points_toric_variety>`
     949          if the default embedding morphism was defined for ``self``,
     950          otherwise a ``ValueError`` exception is raised.
     951
     952        EXAMPLES::
     953
     954            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     955            sage: P1xP1 = ToricVariety(fan, "x s y t")
     956            sage: P1xP1.embedding_morphism()
     957            Traceback (most recent call last):
     958            ...
     959            ValueError: no default embedding was
     960            defined for this toric variety!
     961            sage: patch = P1xP1.affine_patch(0)
     962            sage: patch
     963            2-d affine toric variety
     964            sage: patch.embedding_morphism()
     965            Scheme morphism:
     966              From: 2-d affine toric variety
     967              To:   2-d toric variety covered by 4 affine patches
     968              Defn: Defined on coordinates by sending [x : t] to
     969                    [x : 1 : 1 : t]
     970        """
     971        try:
     972            return self._embedding_morphism
     973        except AttributeError:
     974            raise ValueError("no default embedding was defined for this "
     975                             "toric variety!")
     976
     977    def fan(self, dim=None, codim=None):
     978        r"""
     979        Return the underlying fan of ``self`` or its cones.
     980
     981        INPUT:
     982
     983        - ``dim`` -- dimension of the requested cones;
     984
     985        - ``codim`` -- codimension of the requested cones.
     986
     987        OUTPUT:
     988
     989        - :class:`rational polyhedral fan
     990          <sage.geometry.fan.RationalPolyhedralFan>` if no parameters were
     991          given, :class:`tuple` of :class:`cones
     992          <sage.geometry.cone.ConvexRationalPolyhedralCone>` otherwise.
     993
     994        EXAMPLES::
     995
     996            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     997            sage: P1xP1 = ToricVariety(fan)
     998            sage: P1xP1.fan()
     999            Rational polyhedral fan in 2-d lattice N
     1000            sage: P1xP1.fan() is fan
     1001            True
     1002            sage: P1xP1.fan(1)[0]
     1003            1-d cone of Rational polyhedral fan in 2-d lattice N
     1004        """
     1005        return self._fan(dim, codim)
     1006
     1007    def inject_coefficients(self, scope=None, verbose=True):
     1008        r"""
     1009        Inject generators of the base field of ``self`` into ``scope``.
     1010
     1011        This function is useful if the base field is the field of rational
     1012        functions.
     1013
     1014        INPUT:
     1015
     1016        - ``scope`` -- namespace (default: global);
     1017
     1018        - ``verbose`` -- if ``True`` (default), names of injected generators
     1019          will be printed.
     1020
     1021        OUTPUT:
     1022
     1023        - none.
     1024
     1025        EXAMPLES::
     1026
     1027            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1028            sage: P1xP1 = ToricVariety(fan)
     1029            sage: P1xP1.inject_coefficients()
     1030
     1031        The last command does nothing, since ``P1xP1`` is defined over `\QQ`.
     1032        Let's construct a toric variety over a more complicated field::
     1033
     1034            sage: F = QQ["a, b"].fraction_field()
     1035            sage: P1xP1 = ToricVariety(fan, base_field=F)
     1036            sage: P1xP1.inject_coefficients()
     1037            Defining a, b
     1038        """
     1039        if is_FractionField(self.base_ring()):
     1040            self.base_ring().inject_variables(scope, verbose)
     1041
     1042    def is_homogeneous(self, polynomial):
     1043        r"""
     1044        Check if ``polynomial`` is homogeneous.
     1045
     1046        The coordinate ring of a toric variety is multigraded by relations
     1047        between generating rays of the underlying fan.
     1048
     1049        INPUT:
     1050
     1051        - ``polynomial`` -- polynomial in the coordinate ring of ``self`` or
     1052          its quotient.
     1053
     1054        OUTPUT:
     1055
     1056        - ``True`` if ``polynomial`` is homogeneous and ``False`` otherwise.
     1057
     1058        EXAMPLES:
     1059
     1060        We will use the product of two projective lines with coordinates
     1061        `(x, y)` for one and `(s, t)` for the other::
     1062
     1063            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1064            sage: fan.ray_matrix()
     1065            [ 1  0 -1  0]
     1066            [ 0  1  0 -1]
     1067            sage: P1xP1 = ToricVariety(fan, "x s y t")
     1068            sage: P1xP1.inject_variables()
     1069            Defining x, s, y, t
     1070            sage: P1xP1.is_homogeneous(x - y)
     1071            True
     1072            sage: P1xP1.is_homogeneous(x*s + y*t)
     1073            True
     1074            sage: P1xP1.is_homogeneous(x - t)
     1075            False
     1076        """
     1077        if "_relation_matrix" not in self.__dict__:
     1078            m = self.fan().ray_matrix().transpose().kernel().matrix()
     1079            # We ignore degrees of torus factor coordinates
     1080            m = m.augment(matrix(m.nrows(), self._torus_factor_dim))
     1081            self._relation_matrix = m
     1082        relation_matrix = self._relation_matrix
     1083        if polynomial not in self.coordinate_ring():
     1084            # Then it should be in the quotient corresponding to a subscheme
     1085            polynomial = polynomial.lift()
     1086        monomials = polynomial.monomials()
     1087        if not monomials:
     1088            return True
     1089        degree = relation_matrix * vector(monomials[0].degrees())
     1090        for monomial in monomials:
     1091            if relation_matrix * vector(monomial.degrees()) != degree:
     1092                return False
     1093        return True
     1094
     1095    def is_isomorphic(self, another):
     1096        r"""
     1097        Check if ``self`` is isomorphic to ``another``.
     1098
     1099        INPUT:
     1100
     1101        - ``another`` - :class:`toric variety <ToricVariety_field>`.
     1102
     1103        OUTPUT:
     1104
     1105        - ``True`` if ``self`` and ``another`` are isomorphic,
     1106          ``False`` otherwise.
     1107
     1108        EXAMPLES::
     1109
     1110            sage: fan1 = FaceFan(lattice_polytope.octahedron(2))
     1111            sage: TV1 = ToricVariety(fan1)
     1112            sage: fan2 = NormalFan(lattice_polytope.octahedron(2))
     1113            sage: TV2 = ToricVariety(fan2)
     1114
     1115        Only the most trivial case is implemented so far::
     1116
     1117            sage: TV1.is_isomorphic(TV1)
     1118            True
     1119            sage: TV1.is_isomorphic(TV2)
     1120            Traceback (most recent call last):
     1121            ...
     1122            NotImplementedError:
     1123            isomorphism check is not yet implemented!
     1124        """
     1125        if self is another:
     1126            return True
     1127        if not is_ToricVariety(another):
     1128            raise TypeError(
     1129                "only another toric variety can be checked for isomorphism! "
     1130                "Got %s" % another)
     1131        raise NotImplementedError("isomorphism check is not yet implemented!")
     1132
     1133    def is_complete(self):
     1134        r"""
     1135        Check if ``self`` is complete.
     1136
     1137        OUTPUT:
     1138
     1139        - ``True`` if ``self`` is complete and ``False`` otherwise.
     1140
     1141        EXAMPLES::
     1142
     1143            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1144            sage: P1xP1 = ToricVariety(fan)
     1145            sage: P1xP1.is_complete()
     1146            True
     1147            sage: P1xP1.affine_patch(0).is_complete()
     1148            False
     1149        """
     1150        return self.fan().is_complete()
     1151
     1152    def is_orbifold(self):
     1153        r"""
     1154        Check if ``self`` has only quotient singularities.
     1155
     1156        OUTPUT:
     1157
     1158        - ``True`` if ``self`` has only quotient singularities, ``False``
     1159          otherwise.
     1160
     1161        EXAMPLES::
     1162
     1163            sage: fan1 = FaceFan(lattice_polytope.octahedron(2))
     1164            sage: P1xP1 = ToricVariety(fan1)
     1165            sage: P1xP1.is_orbifold()
     1166            True
     1167            sage: fan2 = NormalFan(lattice_polytope.octahedron(3))
     1168            sage: TV = ToricVariety(fan2)
     1169            sage: TV.is_orbifold()
     1170            False
     1171        """
     1172        return self.fan().is_simplicial()
     1173
     1174    def is_smooth(self):
     1175        r"""
     1176        Check if ``self`` is smooth.
     1177
     1178        OUTPUT:
     1179
     1180        - ``True`` if ``self`` is smooth and ``False`` otherwise.
     1181
     1182        EXAMPLES::
     1183
     1184            sage: diamond = lattice_polytope.octahedron(2)
     1185            sage: fan1 = FaceFan(diamond)
     1186            sage: P1xP1 = ToricVariety(fan1)
     1187            sage: P1xP1.is_smooth()
     1188            True
     1189            sage: fan2 = NormalFan(lattice_polytope.octahedron(2))
     1190            sage: TV = ToricVariety(fan2)
     1191            sage: TV.is_smooth()
     1192            False
     1193        """
     1194        return self.fan().is_smooth()
     1195
     1196    def resolve(self, **kwds):
     1197        r"""
     1198        Construct a toric variety whose fan subdivides the fan of ``self``.
     1199
     1200        The name of this function reflects the fact that usually such
     1201        subdivisions are done for resolving singularities of the original
     1202        variety.
     1203
     1204        INPUT:
     1205
     1206        This function accepts only keyword arguments, none of which are
     1207        mandatory.
     1208
     1209        - ``coordinate_names`` -- names for coordinates of the new variety. If
     1210          not given, will be constructed from the coordinate names of ``self``
     1211          and necessary indexed ones. See :func:`normalize_names` for the
     1212          description of acceptable formats;
     1213
     1214        - ``coordinate_indices`` -- coordinate indices which should be used
     1215          for indexed variables of the new variety;
     1216
     1217        - all other arguments will be passed to
     1218          :meth:`~sage.geometry.fan.RationalPolyhedralFan.subdivide` method of
     1219          the underlying :class:`rational polyhedral fan
     1220          <sage.geometry.fan.RationalPolyhedralFan>`, see its documentation
     1221          for the available options.
     1222
     1223        OUTPUT:
     1224
     1225        - :class:`toric variety <ToricVariety_field>`.
     1226
     1227        EXAMPLES:
     1228
     1229        First we will "manually" resolve a simple orbifold singularity::
     1230
     1231            sage: cone = Cone([(1,1), (-1,1)])
     1232            sage: fan = Fan([cone])
     1233            sage: TV = ToricVariety(fan)
     1234            sage: TV.is_smooth()
     1235            False
     1236            sage: TV_res = TV.resolve(new_rays=[(0,1)])
     1237            sage: TV_res.is_smooth()
     1238            True
     1239            sage: TV_res.fan().ray_matrix()
     1240            [-1  1  0]
     1241            [ 1  1  1]
     1242            sage: [cone.ambient_ray_indices() for cone in TV_res.fan()]
     1243            [(1, 2), (0, 2)]
     1244
     1245        Now let's "automatically" partially resolve a more complicated fan::
     1246
     1247            sage: fan = NormalFan(lattice_polytope.octahedron(3))
     1248            sage: TV = ToricVariety(fan)
     1249            sage: TV.is_smooth()
     1250            False
     1251            sage: TV.is_orbifold()
     1252            False
     1253            sage: TV.fan().nrays()
     1254            8
     1255            sage: TV.fan().ngenerating_cones()
     1256            6
     1257            sage: TV_res = TV.resolve(make_simplicial=True)
     1258            sage: TV_res.is_smooth()
     1259            False
     1260            sage: TV_res.is_orbifold()
     1261            True
     1262            sage: TV_res.fan().nrays()
     1263            8
     1264            sage: TV_res.fan().ngenerating_cones()
     1265            12
     1266            sage: TV.gens()
     1267            (z0, z1, z2, z3, z4, z5, z6, z7)
     1268            sage: TV_res.gens()
     1269            (z0, z1, z2, z3, z4, z5, z6, z7)
     1270            sage: TV_res = TV.resolve(coordinate_names="x+",
     1271            ...                       make_simplicial=True)
     1272            sage: TV_res.gens()
     1273            (x0, x1, x2, x3, x4, x5, x6, x7)
     1274        """
     1275        # If you are changing this function, check out resolve in Fano toric
     1276        # varieties to see if it should be changed too
     1277        #
     1278        # Currently the resolution of fans works for full-dimensional ones
     1279        # only, so there is no point to deal with the general case here, since
     1280        # we will not be able to check that it works.
     1281        coordinate_names = kwds.pop("coordinate_names", None)
     1282        coordinate_indices = kwds.pop("coordinate_indices", None)
     1283        fan = self.fan()
     1284        if fan.dim() != fan.lattice_dim():
     1285            raise NotImplementedError("resolution of toric varieties with "
     1286                                      "torus factors is not yet implemented!")
     1287            # When it is implemented, should be careful with the torus factor
     1288        rfan = fan.subdivide(**kwds)
     1289        if coordinate_names is None:
     1290            coordinate_names = list(self.variable_names())
     1291            if coordinate_indices is None:
     1292                coordinate_indices = range(fan.nrays(), rfan.nrays())
     1293            else:
     1294                coordinate_indices = coordinate_indices[fan.nrays():]
     1295            coordinate_names.extend(normalize_names(
     1296                                    ngens=rfan.nrays() - fan.nrays(),
     1297                                    indices=coordinate_indices,
     1298                                    prefix=self._coordinate_prefix))
     1299            coordinate_names.append(self._coordinate_prefix + "+")
     1300        resolution = ToricVariety(rfan, coordinate_names=coordinate_names,
     1301                                  coordinate_indices=coordinate_indices,
     1302                                  base_field=self.base_ring())
     1303        R = self.coordinate_ring()
     1304        R_res = resolution.coordinate_ring()
     1305        resolution_map = resolution.hom(R.hom(R_res.gens()[:R.ngens()]), self)
     1306        resolution._resolution_map = resolution_map
     1307        # The above map does not have (yet) public methods to access it.
     1308        # While this map is defined correctly, base classes of schemes and
     1309        # morphisms do not treat it as they should. The plan is to fix this
     1310        # situation soon and to be able to use this map!
     1311        return resolution
     1312
     1313    def resolve_to_orbifold(self, **kwds):
     1314        r"""
     1315        Construct an orbifold whose fan subdivides the fan of ``self``.
     1316
     1317        It is a synonym for :meth:`resolve` with ``make_simplicial=True``
     1318        option.
     1319
     1320        INPUT:
     1321
     1322        - this function accepts only keyword arguments. See :meth:`resolve`
     1323          for documentation.
     1324
     1325        OUTPUT:
     1326
     1327        - :class:`toric variety <ToricVariety_field>`.
     1328
     1329        EXAMPLES::
     1330
     1331            sage: fan = NormalFan(lattice_polytope.octahedron(3))
     1332            sage: TV = ToricVariety(fan)
     1333            sage: TV.is_orbifold()
     1334            False
     1335            sage: TV.fan().nrays()
     1336            8
     1337            sage: TV.fan().ngenerating_cones()
     1338            6
     1339            sage: TV_res = TV.resolve_to_orbifold()
     1340            sage: TV_res.is_orbifold()
     1341            True
     1342            sage: TV_res.fan().nrays()
     1343            8
     1344            sage: TV_res.fan().ngenerating_cones()
     1345            12
     1346        """
     1347        return self.resolve(make_simplicial=True, **kwds)
     1348
     1349    def subscheme(self, polynomials):
     1350        r"""
     1351        Return the subscheme of ``self`` defined by ``polynomials``.
     1352
     1353        INPUT:
     1354
     1355        - ``polynomials`` -- list of polynomials in the coordinate ring of
     1356          ``self``.
     1357
     1358        OUTPUT:
     1359
     1360        - :class:`subscheme of a toric variety
     1361          <AlgebraicScheme_subscheme_toric>`.
     1362
     1363        EXAMPLES:
     1364
     1365        We will construct a subscheme of the product of two projective lines
     1366        with coordinates `(x, y)` for one and `(s, t)` for the other::
     1367
     1368            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1369            sage: fan.ray_matrix()
     1370            [ 1  0 -1  0]
     1371            [ 0  1  0 -1]
     1372            sage: P1xP1 = ToricVariety(fan, "x s y t")
     1373            sage: P1xP1.inject_variables()
     1374            Defining x, s, y, t
     1375            sage: X = P1xP1.subscheme([x*s + y*t, x^3+y^3])
     1376            sage: X
     1377            Closed subscheme of 2-d toric variety
     1378            covered by 4 affine patches defined by:
     1379              x*s + y*t,
     1380              x^3 + y^3
     1381            sage: X.defining_polynomials()
     1382            (x*s + y*t, x^3 + y^3)
     1383            sage: X.defining_ideal()
     1384            Ideal (x*s + y*t, x^3 + y^3)
     1385            of Multivariate Polynomial Ring in x, s, y, t
     1386            over Rational Field
     1387            sage: X.base_ring()
     1388            Rational Field
     1389            sage: X.base_scheme()
     1390            Spectrum of Rational Field
     1391            sage: X.structure_morphism()
     1392            Scheme morphism:
     1393              From: Closed subscheme of
     1394              2-d toric variety covered by 4 affine patches defined by:
     1395              x*s + y*t,
     1396              x^3 + y^3
     1397              To:   Spectrum of Rational Field
     1398              Defn: Structure map
     1399        """
     1400        return AlgebraicScheme_subscheme_toric(self, polynomials)
     1401
     1402
     1403class AlgebraicScheme_subscheme_toric(AlgebraicScheme_subscheme):
     1404    r"""
     1405    Construct an algebraic subscheme of a toric variety.
     1406
     1407    .. WARNING::
     1408
     1409        You should not create objects of this class directly. The preferred
     1410        method to construct such subschemes is to use
     1411        :meth:`~ToricVariety_field.subscheme` method of
     1412        :class:`toric varieties <ToricVariety_field>`.
     1413
     1414    INPUT:
     1415
     1416    - ``toric_variety`` -- ambient :class:`toric variety
     1417      <ToricVariety_field>`;
     1418
     1419    - ``polynomials`` -- single polynomial, list, or ideal of defining
     1420      polynomials in the coordinate ring of ``toric_variety``.
     1421
     1422    OUTPUT:
     1423
     1424    - :class:`algebraic subscheme of a toric variety
     1425      <AlgebraicScheme_subscheme_toric>`.
     1426
     1427    TESTS::
     1428
     1429        sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1430        sage: P1xP1 = ToricVariety(fan, "x s y t")
     1431        sage: P1xP1.inject_variables()
     1432        Defining x, s, y, t
     1433        sage: import sage.schemes.generic.toric_variety as tv
     1434        sage: X = tv.AlgebraicScheme_subscheme_toric(
     1435        ...         P1xP1, [x*s + y*t, x^3+y^3])
     1436        sage: X
     1437        Closed subscheme of 2-d toric variety
     1438        covered by 4 affine patches defined by:
     1439          x*s + y*t,
     1440          x^3 + y^3
     1441
     1442    A better way to construct the same scheme as above::
     1443
     1444        sage: P1xP1.subscheme([x*s + y*t, x^3+y^3])
     1445        Closed subscheme of 2-d toric variety
     1446        covered by 4 affine patches defined by:
     1447          x*s + y*t,
     1448          x^3 + y^3
     1449    """
     1450
     1451    def __init__(self, toric_variety, polynomials):
     1452        r"""
     1453        See :class:`AlgebraicScheme_subscheme_toric` for documentation.
     1454
     1455        TESTS::
     1456
     1457            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1458            sage: P1xP1 = ToricVariety(fan, "x s y t")
     1459            sage: P1xP1.inject_variables()
     1460            Defining x, s, y, t
     1461            sage: import sage.schemes.generic.toric_variety as tv
     1462            sage: X = tv.AlgebraicScheme_subscheme_toric(
     1463            ...         P1xP1, [x*s + y*t, x^3+y^3])
     1464            sage: X
     1465            Closed subscheme of 2-d toric variety
     1466            covered by 4 affine patches defined by:
     1467              x*s + y*t,
     1468              x^3 + y^3
     1469        """
     1470        # Just to make sure that keyword arguments will be passed correctly
     1471        super(AlgebraicScheme_subscheme_toric, self).__init__(toric_variety,
     1472                                                              polynomials)
     1473
     1474    def _point_morphism_class(self, *args, **kwds):
     1475        r"""
     1476        Construct a morphism determined by action on points of ``self``.
     1477
     1478        INPUT:
     1479
     1480        - same as for :class:`SchemeMorphism_on_points_toric_variety`.
     1481
     1482        OUPUT:
     1483
     1484        - :class:`SchemeMorphism_on_points_toric_variety`.
     1485
     1486        TESTS::
     1487
     1488            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1489            sage: P1xP1 = ToricVariety(fan)
     1490            sage: P1xP1.inject_variables()
     1491            Defining z0, z1, z2, z3
     1492            sage: P1 = P1xP1.subscheme(z0-z2)
     1493            sage: H = P1.Hom(P1xP1)
     1494            sage: P1._point_morphism_class(H, [z0,z1,z0,z3])
     1495            Scheme morphism:
     1496              From: Closed subscheme of 2-d toric variety
     1497              covered by 4 affine patches defined by:
     1498              z0 - z2
     1499              To:   2-d toric variety covered by 4 affine patches
     1500              Defn: Defined on coordinates by sending [z0 : z1 : z2 : z3] to
     1501                    [z2bar : z1bar : z2bar : z3bar]
     1502        """
     1503        return SchemeMorphism_on_points_toric_variety(*args, **kwds)
     1504
     1505    def affine_patch(self, i):
     1506        r"""
     1507        Return the ``i``-th affine patch of ``self``.
     1508
     1509        INPUT:
     1510
     1511        - ``i`` -- integer, index of a generating cone of the fan of the
     1512          ambient space of ``self``.
     1513
     1514        OUTPUT:
     1515
     1516        - subscheme of an affine :class:`toric variety <ToricVariety_field>`
     1517          corresponding to the pull-back of ``self`` by the embedding morphism
     1518          of the ``i``-th :meth:`affine patch of the ambient space
     1519          <ToricVariety_field.affine_patch>` of ``self``.
     1520
     1521        The result is cached, so the ``i``-th patch is always the same object
     1522        in memory.
     1523
     1524        EXAMPLES::
     1525
     1526            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1527            sage: P1xP1 = ToricVariety(fan, "x s y t")
     1528            sage: patch1 = P1xP1.affine_patch(1)
     1529            sage: patch1.embedding_morphism()
     1530            Scheme morphism:
     1531              From: 2-d affine toric variety
     1532              To:   2-d toric variety covered by 4 affine patches
     1533              Defn: Defined on coordinates by sending [y : t] to
     1534                    [1 : 1 : y : t]
     1535            sage: P1xP1.inject_variables()
     1536            Defining x, s, y, t
     1537            sage: P1 = P1xP1.subscheme(x-y)
     1538            sage: subpatch = P1.affine_patch(1)
     1539            sage: subpatch
     1540            Closed subscheme of 2-d affine toric variety defined by:
     1541              -y + 1
     1542        """
     1543        i = int(i)   # implicit type checking
     1544        try:
     1545            return self._affine_patches[i]
     1546        except AttributeError:
     1547            self._affine_patches = dict()
     1548        except KeyError:
     1549            pass
     1550        ambient_patch = self.ambient_space().affine_patch(i)
     1551        phi_p = ambient_patch.embedding_morphism().defining_polynomials()
     1552        patch = ambient_patch.subscheme(
     1553                            [p(phi_p) for p in self.defining_polynomials()])
     1554        patch._embedding_morphism = patch.hom(phi_p, self)
     1555        self._affine_patches[i] = patch
     1556        return patch
     1557
     1558    def dimension(self):
     1559        """
     1560        Return the dimension of ``self``.
     1561
     1562        .. NOTE::
     1563
     1564            Currently the dimension of subschemes of toric varieties can be
     1565            returned only if it was somehow set before.
     1566
     1567        OUTPUT:
     1568
     1569        - integer.
     1570
     1571        EXAMPLES::
     1572
     1573            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1574            sage: P1xP1 = ToricVariety(fan)
     1575            sage: P1xP1.inject_variables()
     1576            Defining z0, z1, z2, z3
     1577            sage: P1 = P1xP1.subscheme(z0-z2)
     1578            sage: P1.dimension()
     1579            Traceback (most recent call last):
     1580            ...
     1581            NotImplementedError:
     1582            cannot compute dimension of this scheme!
     1583        """
     1584        if "_dimension" not in self.__dict__:
     1585            raise NotImplementedError(
     1586                                "cannot compute dimension of this scheme!")
     1587        return self._dimension
     1588
     1589    def embedding_morphism(self):
     1590        r"""
     1591        Return the default embedding morphism of ``self``.
     1592
     1593        Such a morphism is always defined for an affine patch of a subscheme
     1594        of a toric variety (which is a subscheme of a toric variety itself).
     1595
     1596        OUTPUT:
     1597
     1598        - :class:`scheme morphism <SchemeMorphism_on_points_toric_variety>`
     1599          if the default embedding morphism was defined for ``self``,
     1600          otherwise a ``ValueError`` exception is raised.
     1601
     1602        EXAMPLES::
     1603
     1604            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1605            sage: P1xP1 = ToricVariety(fan, "x s y t")
     1606            sage: patch1 = P1xP1.affine_patch(1)
     1607            sage: patch1.embedding_morphism()
     1608            Scheme morphism:
     1609              From: 2-d affine toric variety
     1610              To:   2-d toric variety covered by 4 affine patches
     1611              Defn: Defined on coordinates by sending [y : t] to
     1612                    [1 : 1 : y : t]
     1613            sage: P1xP1.inject_variables()
     1614            Defining x, s, y, t
     1615            sage: P1 = P1xP1.subscheme(x-y)
     1616            sage: P1.embedding_morphism()
     1617            Traceback (most recent call last):
     1618            ...
     1619            ValueError: no default embedding was defined
     1620            for this subscheme of a toric variety!
     1621            sage: subpatch = P1.affine_patch(1)
     1622            sage: subpatch
     1623            Closed subscheme of 2-d affine toric variety defined by:
     1624              -y + 1
     1625            sage: subpatch.embedding_morphism()
     1626            Scheme morphism:
     1627              From: Closed subscheme of 2-d affine toric variety defined by:
     1628              -y + 1
     1629              To:   Closed subscheme of 2-d toric variety
     1630              covered by 4 affine patches defined by:
     1631              x - y
     1632              Defn: Defined on coordinates by sending [y : t] to
     1633                    [1 : 1 : 1 : tbar]
     1634        """
     1635        try:
     1636            return self._embedding_morphism
     1637        except AttributeError:
     1638            raise ValueError("no default embedding was defined for this "
     1639                             "subscheme of a toric variety!")
     1640
     1641
     1642class SchemeHomset_toric_coordinates_field(SchemeHomset_coordinates):
     1643    """
     1644    Construct the `Hom`-space of morphisms given on coordinates.
     1645
     1646    .. WARNING::
     1647
     1648        You should not create objects of this class directly.
     1649
     1650
     1651    INPUT:
     1652
     1653    - same as for :class:`sage.schemes.generic.SchemeHomset_coordinates`.
     1654
     1655    OUPUT:
     1656
     1657    - :class:`SchemeHomset_toric_coordinates_field`.
     1658
     1659    TESTS::
     1660
     1661        sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1662        sage: P1xP1 = ToricVariety(fan)
     1663        sage: import sage.schemes.generic.toric_variety as tv
     1664        sage: tv.SchemeHomset_toric_coordinates_field(P1xP1, QQ)
     1665        Set of Rational Points of 2-d toric variety
     1666        covered by 4 affine patches
     1667
     1668    A better way to construct the same `Hom`-space as above::
     1669
     1670        sage: P1xP1(QQ)
     1671        Set of Rational Points of 2-d toric variety
     1672        covered by 4 affine patches
     1673    """
     1674    # Mimicking SchemeHomset_projective_coordinates_field,
     1675    # affine spaces implement only "except" case
     1676    def __call__(self, *arg):
     1677        r"""
     1678        Construct a morphism from given parameters.
     1679
     1680        INPUT:
     1681
     1682        - data determining a morphism.
     1683
     1684        TESTS::
     1685
     1686            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1687            sage: P1xP1 = ToricVariety(fan)
     1688            sage: import sage.schemes.generic.toric_variety as tv
     1689            sage: H = tv.SchemeHomset_toric_coordinates_field(P1xP1, QQ)
     1690            sage: H(1,2,3,4)
     1691            [1 : 2 : 3 : 4]
     1692        """
     1693        # This may break for one-dimensional varieties.
     1694        if len(arg) == 1:
     1695            arg = arg[0]
     1696        X = self.codomain()
     1697        try:
     1698            return X._point_class(X, arg)
     1699        except AttributeError:  # should be very rare
     1700            return SchemeMorphism_toric_coordinates_field(self, arg)
     1701
     1702
     1703class SchemeMorphism_toric_coordinates_field(SchemeMorphism_coordinates):
     1704    """
     1705    Construct a morphism determined by giving coordinates in a field.
     1706
     1707    .. WARNING::
     1708
     1709        You should not create objects of this class directly.
     1710
     1711    INPUT:
     1712
     1713    - ``X`` -- subscheme of a toric variety.
     1714
     1715    - ``coordinates`` -- list of coordinates in the base field of ``X``.
     1716
     1717    - ``check`` -- if ``True`` (default), the input will be checked for
     1718      correctness.
     1719
     1720    OUTPUT:
     1721
     1722    - :class:`SchemeMorphism_toric_coordinates_field`.
     1723
     1724    TESTS::
     1725
     1726        sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1727        sage: P1xP1 = ToricVariety(fan)
     1728        sage: P1xP1(1,2,3,4)
     1729        [1 : 2 : 3 : 4]
     1730    """
     1731    # Mimicking affine/projective classes
     1732    def __init__(self, X, coordinates, check=True):
     1733        r"""
     1734        See :class:`SchemeMorphism_toric_coordinates_field` for documentation.
     1735
     1736        TESTS::
     1737
     1738            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1739            sage: P1xP1 = ToricVariety(fan)
     1740            sage: P1xP1(1,2,3,4)
     1741            [1 : 2 : 3 : 4]
     1742        """
     1743        # Convert scheme to its set of points over the base ring
     1744        if is_Scheme(X):
     1745            X = X(X.base_ring())
     1746        super(SchemeMorphism_toric_coordinates_field, self).__init__(X)
     1747        if check:
     1748            # Verify that there are the right number of coords
     1749            # Why is it not done in the parent?
     1750            if is_SchemeMorphism(coordinates):
     1751                coordinates = list(coordinates)
     1752            if not isinstance(coordinates, (list, tuple)):
     1753                raise TypeError("coordinates must be a scheme point, list, "
     1754                                "or tuple. Got %s" % coordinates)
     1755            d = X.codomain().ambient_space().ngens()
     1756            if len(coordinates) != d:
     1757                raise ValueError("there must be %d coordinates! Got only %d: "
     1758                                 "%s" % (d, len(coordinates), coordinates))
     1759            # Make sure the coordinates all lie in the appropriate ring
     1760            coordinates = Sequence(coordinates, X.value_ring())
     1761            # Verify that the point satisfies the equations of X.
     1762            X.codomain()._check_satisfies_equations(coordinates)
     1763        self._coords = coordinates
     1764
     1765
     1766class SchemeMorphism_on_points_toric_variety(SchemeMorphism_on_points):
     1767    """
     1768    Construct a morphism determined by action on points.
     1769
     1770    .. WARNING::
     1771
     1772        You should not create objects of this class directly.
     1773
     1774    INPUT:
     1775
     1776    - same as for :class:`sage.schemes.generic.SchemeMorphism_on_points`.
     1777
     1778    OUPUT:
     1779
     1780    - :class:`SchemeMorphism_on_points_toric_variety`.
     1781
     1782    TESTS::
     1783
     1784        sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1785        sage: P1xP1 = ToricVariety(fan)
     1786        sage: P1xP1.inject_variables()
     1787        Defining z0, z1, z2, z3
     1788        sage: P1 = P1xP1.subscheme(z0-z2)
     1789        sage: H = P1xP1.Hom(P1)
     1790        sage: import sage.schemes.generic.toric_variety as tv
     1791        sage: tv.SchemeMorphism_on_points_toric_variety(H, [z0,z1,z0,z3])
     1792        Scheme morphism:
     1793          From: 2-d toric variety covered by 4 affine patches
     1794          To:   Closed subscheme of 2-d toric variety
     1795                covered by 4 affine patches defined by:
     1796          z0 - z2
     1797          Defn: Defined on coordinates by sending
     1798                [z0 : z1 : z2 : z3] to [z0 : z1 : z0 : z3]
     1799    """
     1800
     1801    def __init__(self, parent, polynomials, check=True):
     1802        r"""
     1803        See :class:`SchemeMorphism_on_points_toric_variety` for documentation.
     1804
     1805        TESTS::
     1806
     1807            sage: fan = FaceFan(lattice_polytope.octahedron(2))
     1808            sage: P1xP1 = ToricVariety(fan)
     1809            sage: P1xP1.inject_variables()
     1810            Defining z0, z1, z2, z3
     1811            sage: P1 = P1xP1.subscheme(z0-z2)
     1812            sage: H = P1xP1.Hom(P1)
     1813            sage: import sage.schemes.generic.toric_variety as tv
     1814            sage: tv.SchemeMorphism_on_points_toric_variety(H, [z0,z1,z0,z3])
     1815            Scheme morphism:
     1816              From: 2-d toric variety covered by 4 affine patches
     1817              To:   Closed subscheme of 2-d toric variety
     1818                    covered by 4 affine patches defined by:
     1819              z0 - z2
     1820              Defn: Defined on coordinates by sending
     1821                    [z0 : z1 : z2 : z3] to [z0 : z1 : z0 : z3]
     1822        """
     1823        SchemeMorphism_on_points.__init__(self, parent, polynomials, check)
     1824        if check:
     1825            # Check that defining polynomials are homogeneous (degrees can be
     1826            # different if the target uses weighted coordinates)
     1827            for p in self.defining_polynomials():
     1828                if not self.domain().ambient_space().is_homogeneous(p):
     1829                    raise ValueError("%s is not homogeneous!" % p)
     1830
     1831
     1832def normalize_names(names=None, ngens=None, prefix=None, indices=None,
     1833                    return_prefix=False):
     1834    r"""
     1835    Return a list of names in the standard form.
     1836
     1837    INPUT:
     1838
     1839    All input parameters are optional.
     1840
     1841    - ``names`` -- names given either as a single string (with individual
     1842      names separated by commas or spaces) or a list of strings with each
     1843      string specifying a name. If the last name ends with the plus sign,
     1844      "+", this name will be used as ``prefix`` (even if ``prefix`` was
     1845      given explicitly);
     1846
     1847    - ``ngens`` -- number of names to be returned;
     1848
     1849    - ``prefix`` -- prefix for the indexed names given as a string;
     1850
     1851    - ``indices`` -- list of integers (default: ``range(ngens)``) used as
     1852      indices for names with ``prefix``. If given, must be of length
     1853      ``ngens``;
     1854
     1855    - ``return_prefix`` -- if ``True``, the last element of the returned list
     1856      will contain the prefix determined from ``names`` or given as the
     1857      parameter ``prefix``. This is useful if you may need more names in the
     1858      future.
     1859
     1860    OUTPUT:
     1861
     1862    - list of names given as strings.
     1863
     1864    These names are constructed in the following way:
     1865
     1866    #. If necessary, split ``names`` into separate names.
     1867    #. If the last name ends with "+", put it into ``prefix``.
     1868    #. If ``ngens`` was given, add to the names obtained so far as many
     1869       indexed names as necessary to get this number. If the ``k``-th name of
     1870       the *total* list of names is indexed, it is
     1871       ``prefix + str(indices[k])``. If there were already more names than
     1872       ``ngens``, discard "extra" ones.
     1873    #. Check if constructed names are valid. See :func:`certify_names` for
     1874       details.
     1875    #. If the option ``return_prefix=True`` was given, add ``prefix`` to the
     1876       end of the list.
     1877
     1878    EXAMPLES:
     1879
     1880    As promised, all parameters are optional::
     1881
     1882        sage: from sage.schemes.generic.toric_variety import normalize_names
     1883        sage: normalize_names()
     1884        []
     1885
     1886    One of the most common uses is probably this one::
     1887
     1888        sage: normalize_names("x+", 4)
     1889        ['x0', 'x1', 'x2', 'x3']
     1890
     1891    Now suppose that you want to enumerate your variables starting with one
     1892    instead of zero::
     1893
     1894        sage: normalize_names("x+", 4, indices=range(1,5))
     1895        ['x1', 'x2', 'x3', 'x4']
     1896
     1897    You may actually have an arbitrary enumeration scheme::
     1898
     1899        sage: normalize_names("x+", 4, indices=[1, 10, 100, 1000])
     1900        ['x1', 'x10', 'x100', 'x1000']
     1901
     1902    Now let's add some "explicit" names::
     1903
     1904        sage: normalize_names("x y z t+", 4)
     1905        ['x', 'y', 'z', 't3']
     1906
     1907    Note that the "automatic" name is ``t3`` instead of ``t0``. This may seem
     1908    weird, but the reason for this behaviour is that the fourth name in this
     1909    list will be the same no matter how many explicit names were given::
     1910
     1911        sage: normalize_names("x y t+", 4)
     1912        ['x', 'y', 't2', 't3']
     1913
     1914    This is especially useful if you get ``names`` from a user but want to
     1915    specify all default names::
     1916
     1917        sage: normalize_names("x, y", 4, prefix="t")
     1918        ['x', 'y', 't2', 't3']
     1919
     1920    In this format, the user can easily override your choice for automatic
     1921    names::
     1922
     1923        sage: normalize_names("x y s+", 4, prefix="t")
     1924        ['x', 'y', 's2', 's3']
     1925
     1926    Let's now use all parameters at once::
     1927
     1928        sage: normalize_names("x, y, s+", 4, prefix="t",
     1929        ...       indices=range(1,5), return_prefix=True)
     1930        ['x', 'y', 's3', 's4', 's']
     1931
     1932    Note that you still need to give indices for all names, even if some of
     1933    the first ones will be "wasted" because of the explicit names. The reason
     1934    is the same as before - this ensures consistency of automatically
     1935    generated names, no matter how many explicit names were given.
     1936
     1937    The prefix is discarded if ``ngens`` was not given::
     1938
     1939        sage: normalize_names("alpha, beta, gamma, zeta+")
     1940        ['alpha', 'beta', 'gamma']
     1941
     1942    Finally, let's take a look at some possible mistakes::
     1943
     1944        sage: normalize_names("123")
     1945        Traceback (most recent call last):
     1946        ...
     1947        ValueError: name must start with a letter! Got 123
     1948
     1949    A more subtle one::
     1950
     1951        sage: normalize_names("x1", 4, prefix="x")
     1952        Traceback (most recent call last):
     1953        ...
     1954        ValueError: names must be distinct! Got: ['x1', 'x1', 'x2', 'x3']
     1955    """
     1956    if names is None:
     1957        names = []
     1958    elif isinstance(names, str):
     1959        names = names.replace(",", " ").split()
     1960    else:
     1961        try:
     1962            names = list(names)
     1963        except TypeError:
     1964            raise TypeError(
     1965                    "names must be a string or a list or tuple of them!")
     1966        for name in names:
     1967            if not isinstance(name, str):
     1968                raise TypeError(
     1969                    "names must be a string or a list or tuple of them!")
     1970    if names and names[-1].endswith("+"):
     1971        prefix = names.pop()[:-1]
     1972    if ngens is None:
     1973        ngens = len(names)
     1974    if len(names) < ngens:
     1975        if prefix is None:
     1976            raise IndexError("need %d names but only %d are given!"
     1977                             % (ngens, len(names)))
     1978        if indices is None:
     1979            indices = range(ngens)
     1980        elif len(indices) != ngens:
     1981            raise ValueError("need exactly %d indices, but got %d!"
     1982                             % (ngens, len(indices)))
     1983        names += [prefix + str(i) for i in indices[len(names):]]
     1984    if len(names) > ngens:
     1985        names = names[:ngens]
     1986    # Check that all given and constructed names are valid
     1987    certify_names(names)
     1988    if return_prefix:
     1989        names.append(prefix)
     1990    return names
     1991
     1992
     1993def certify_names(names):
     1994    r"""
     1995    Make sure that ``names`` are valid in Python.
     1996
     1997    INPUT:
     1998
     1999    - ``names`` -- list of strings.
     2000
     2001    OUTPUT:
     2002
     2003    - none, but a ``ValueError`` exception is raised if ``names`` are invalid.
     2004
     2005    Each name must satisfy the following requirements:
     2006
     2007    * Be non-empty.
     2008    * Contain only (Latin) letters, digits, and underscores ("_").
     2009    * Start with a letter.
     2010
     2011    In addition, all names must be distinct.
     2012
     2013    EXAMPLES::
     2014
     2015        sage: from sage.schemes.generic.toric_variety import certify_names
     2016        sage: certify_names([])
     2017        sage: certify_names(["a", "x0", "x_45"])
     2018        sage: certify_names(["", "x0", "x_45"])
     2019        Traceback (most recent call last):
     2020        ...
     2021        ValueError: name must be nonempty!
     2022        sage: certify_names(["a", "0", "x_45"])
     2023        Traceback (most recent call last):
     2024        ...
     2025        ValueError: name must start with a letter! Got 0
     2026        sage: certify_names(["a", "x0", "Щ_45"])
     2027        Traceback (most recent call last):
     2028        ...
     2029        ValueError: name must be alphanumeric! Got Щ_45
     2030        sage: certify_names(["a", "x0", "x0"])
     2031        Traceback (most recent call last):
     2032        ...
     2033        ValueError: names must be distinct! Got: ['a', 'x0', 'x0']
     2034    """
     2035    for name in names:
     2036        if not name:
     2037            raise ValueError("name must be nonempty!")
     2038        if not name.isalnum() and not name.replace("_","").isalnum():
     2039            # Must be alphanumeric except for non-leading '_'
     2040            raise ValueError("name must be alphanumeric! Got %s" % name)
     2041        if not name[0].isalpha():
     2042            raise ValueError("name must start with a letter! Got %s" % name)
     2043    if len(set(names)) != len(names):
     2044        raise ValueError("names must be distinct! Got: %s" % names)