Ticket #9972: trac_9972_add_fan_morphisms.patch

File trac_9972_add_fan_morphisms.patch, 37.2 KB (added by novoselt, 9 years ago)

preimage_cones implemented

  • doc/en/reference/geometry.rst

    # HG changeset patch
    # User Andrey Novoseltsev <novoselt@gmail.com>
    # Date 1289429460 25200
    # Node ID 98a87b6883296d2bea0c123b81b42c17dcce5809
    # Parent  8d193571aaec39f6c2f924baeba0d5b94a834ddf
    Trac 9972: Add support for fan morphisms.
    
    diff -r 8d193571aaec -r 98a87b688329 doc/en/reference/geometry.rst
    a b  
    88.. toctree::
    99   :maxdepth: 2
    1010
     11   sage/geometry/toric_lattice   
    1112   sage/geometry/cone
     13   sage/geometry/fan
     14   sage/geometry/fan_morphism
     15   sage/geometry/toric_plotter
     16   
    1217   sage/rings/polynomial/groebner_fan
     18   
    1319   sage/geometry/lattice_polytope
    1420   sage/geometry/polyhedra
    15    sage/geometry/fan
    16    
    17    sage/geometry/toric_lattice
    18    
    19    sage/geometry/toric_plotter
    2021
  • sage/geometry/all.py

    diff -r 8d193571aaec -r 98a87b688329 sage/geometry/all.py
    a b  
    22
    33from fan import Fan, FaceFan, NormalFan
    44
     5from fan_morphism import FanMorphism
     6
    57from polytope import polymake
    68
    79from polyhedra import Polyhedron, polytopes
  • new file sage/geometry/fan_morphism.py

    diff -r 8d193571aaec -r 98a87b688329 sage/geometry/fan_morphism.py
    - +  
     1r"""
     2Morphisms between toric lattices compatible with fans
     3
     4This module is a part of the framework for toric varieties
     5(:mod:`~sage.schemes.generic.toric_variety`,
     6:mod:`~sage.schemes.generic.fano_toric_variety`). Its main purpose is to
     7provide support for working with lattice morphisms compatible with fans via
     8:class:`FanMorphism` class.
     9
     10AUTHORS:
     11
     12- Andrey Novoseltsev (2010-10-17): initial version.
     13
     14EXAMPLES:
     15
     16Let's consider the face and normal fans of the "diamond" and the projection
     17to the `x`-axis::
     18
     19    sage: diamond = lattice_polytope.octahedron(2)
     20    sage: face = FaceFan(diamond)
     21    sage: normal = NormalFan(diamond)
     22    sage: N = face.lattice()
     23    sage: H = End(N)
     24    sage: phi = H([N.0, 0])
     25    sage: phi     
     26    Free module morphism defined by the matrix
     27    [1 0]
     28    [0 0]
     29    Domain: 2-d lattice N
     30    Codomain: 2-d lattice N
     31    sage: FanMorphism(phi, normal, face)
     32    Traceback (most recent call last):
     33    ...
     34    ValueError: the image of generating cone #1 of the domain fan
     35    is not contained in a single cone of the codomain fan!
     36   
     37Some of the cones of the normal fan fail to be mapped to a single cone of the
     38face fan. We can rectify the situation in the following way::
     39
     40    sage: fm = FanMorphism(phi, normal, face, subdivide=True)
     41    sage: fm
     42    Fan morphism defined by the matrix
     43    [1 0]
     44    [0 0]
     45    Domain fan: Rational polyhedral fan in 2-d lattice N
     46    Codomain fan: Rational polyhedral fan in 2-d lattice N   
     47    sage: fm.domain_fan().ray_matrix()
     48    [-1  1 -1  1  0  0]
     49    [ 1  1 -1 -1 -1  1]   
     50    sage: normal.ray_matrix()
     51    [-1  1 -1  1]
     52    [ 1  1 -1 -1]
     53   
     54As you see, it was necessary to insert two new rays (to prevent "upper" and
     55"lower" cones of the normal fan from being mapped to the whole `x`-axis).
     56"""
     57
     58
     59#*****************************************************************************
     60#       Copyright (C) 2010 Andrey Novoseltsev <novoselt@gmail.com>
     61#       Copyright (C) 2010 William Stein <wstein@gmail.com>
     62#
     63#  Distributed under the terms of the GNU General Public License (GPL)
     64#
     65#                  http://www.gnu.org/licenses/
     66#*****************************************************************************
     67
     68
     69import operator
     70
     71from sage.categories.all import Hom
     72from sage.geometry.cone import Cone
     73from sage.geometry.fan import Fan, is_Fan
     74from sage.matrix.all import is_Matrix
     75from sage.misc.all import latex, walltime
     76from sage.modules.free_module_morphism import (FreeModuleMorphism,
     77                                               is_FreeModuleMorphism)
     78
     79
     80class FanMorphism(FreeModuleMorphism):
     81    r"""
     82    Create a fan morphism.
     83   
     84    Let `\Sigma_1` and `\Sigma_2` be two fans in lattices `N_1` and `N_2`
     85    respectively. Let `\phi` be a morphism (i.e. a linear map) from `N_1` to
     86    `N_2`. We say that `\phi` is *compatible* with `\Sigma_1` and `\Sigma_2`
     87    if every cone `\sigma_1\in\Sigma_1` is mapped by `\phi` into a single cone
     88    `\sigma_2\in\Sigma_2`, i.e. `\phi(\sigma_1)\subset\sigma_2` (`\sigma_2`
     89    may be different for different `\sigma_1`).
     90   
     91    By a **fan morphism** we understand a morphism between two lattices
     92    compatible with specified fans in these lattices. Such morphisms behave in
     93    exactly the same way as "regular" morphisms between lattices, but:
     94   
     95    * fan morphisms have a special constructor allowing some automatic
     96      adjustments to the initial fans (see below);
     97    * fan morphisms are aware of the associated fans and they can be
     98      accessed via :meth:`codomain_fan` and :meth:`domain_fan`;
     99    * fan morphisms can efficiently compute :meth:`image_cone` of a given
     100      cone of the domain fan and :meth:`preimage_cones` of a given cone of
     101      the codomain fan.
     102   
     103    INPUT:
     104       
     105    - ``morphism`` -- either a morphism between domain and codomain, or an
     106      integral matrix defining such a morphism;
     107     
     108    - ``domain_fan`` -- a :class:`fan
     109      <sage.geometry.fan.RationalPolyhedralFan>` in the domain;
     110     
     111    - ``codomain`` -- (default: ``None``) either a codomain lattice or a fan in
     112      the codomain. If the codomain fan is not given, the image fan (fan
     113      generated by images of generating cones) of ``domain_fan`` will be used,
     114      if possible;
     115     
     116    - ``subdivide`` -- (default: ``False``) if ``True`` and ``domain_fan`` is
     117      not compatible with the codomain fan because it is too coarse, it will be
     118      automatically refined to become compatible (the minimal refinement is
     119      canonical, so there are no choices involved);
     120     
     121    - ``check`` -- (default: ``True``) if ``False``, given fans and morphism
     122      will be assumed to be compatible. Be careful when using this option,
     123      since wrong assumptions can lead to wrong and hard-to-detect errors. On
     124      the other hand, this option may save you some time;
     125     
     126    - ``verbose`` -- (default: ``False``) if ``True``, some information may be
     127      printed during construction of the fan morphism.
     128
     129    OUTPUT:
     130   
     131    - a fan morphism.
     132   
     133    EXAMPLES:
     134   
     135    Here we consider the face and normal fans of the "diamond" and the
     136    projection to the `x`-axis::
     137   
     138        sage: diamond = lattice_polytope.octahedron(2)
     139        sage: face = FaceFan(diamond)
     140        sage: normal = NormalFan(diamond)
     141        sage: N = face.lattice()
     142        sage: H = End(N)
     143        sage: phi = H([N.0, 0])
     144        sage: phi     
     145        Free module morphism defined by the matrix
     146        [1 0]
     147        [0 0]
     148        Domain: 2-d lattice N
     149        Codomain: 2-d lattice N
     150        sage: fm = FanMorphism(phi, face, normal)
     151        sage: fm.domain_fan() is face
     152        True
     153       
     154    Note, that since ``phi`` is compatible with these fans, the returned
     155    fan is exactly the same object as the initial ``domain_fan``. ::
     156   
     157        sage: FanMorphism(phi, normal, face)
     158        Traceback (most recent call last):
     159        ...
     160        ValueError: the image of generating cone #1 of the domain fan
     161        is not contained in a single cone of the codomain fan!           
     162        sage: fm = FanMorphism(phi, normal, face, subdivide=True)
     163        sage: fm.domain_fan() is normal
     164        False
     165        sage: fm.domain_fan().ngenerating_cones()
     166        6
     167       
     168    We had to subdivide two of the four cones of the normal fan, since
     169    they were mapped by ``phi`` into non-strictly convex cones.
     170   
     171    It is possible to omit the codomain fan, in which case the image fan will
     172    be used instead of it::
     173   
     174        sage: fm = FanMorphism(phi, face)
     175        sage: fm.codomain_fan()
     176        Rational polyhedral fan in 2-d lattice N
     177        sage: fm.codomain_fan().ray_matrix()
     178        [ 1 -1]
     179        [ 0  0]
     180       
     181    Now we demonstrate a more subtle example. We take the first quadrant as our
     182    domain fan. Then we divide the first quadrant into three cones, throw away
     183    the middle one and take the other two as our codomain fan. These fans are
     184    incompatible with the identity lattice morphism since the image of the
     185    domain fan is out of the support of the codomain fan::
     186   
     187        sage: N = ToricLattice(2)
     188        sage: phi = End(N).identity()
     189        sage: F1 = Fan(cones=[(0,1)], rays=[(1,0), (0,1)])
     190        sage: F2 = Fan(cones=[(0,1), (2,3)],
     191        ...            rays=[(1,0), (2,1), (1,2), (0,1)])
     192        sage: FanMorphism(phi, F1, F2)
     193        Traceback (most recent call last):
     194        ...
     195        ValueError: the image of generating cone #0 of the domain fan
     196        is not contained in a single cone of the codomain fan!
     197        sage: FanMorphism(phi, F1, F2, subdivide=True)
     198        Traceback (most recent call last):
     199        ...
     200        ValueError: morphism defined by
     201        [1 0]
     202        [0 1]
     203        does not map
     204        Rational polyhedral fan in 2-d lattice N
     205        into the support of
     206        Rational polyhedral fan in 2-d lattice N!
     207       
     208    The problem was detected and handled correctly (i.e. an exception was
     209    raised). However, the used algorithm requires extra checks for this
     210    situation after constructing a potential subdivision and this can take
     211    significant time. You can save about half the time using
     212    ``check=False`` option, if you know in advance that it is possible to
     213    make fans compatible with the morphism by subdividing the domain fan.
     214    Of course, if your assumption was incorrect, the result will be wrong
     215    and you will get a fan which *does* map into the support of the
     216    codomain fan, but is **not** a subdivision of the domain fan. You
     217    can test it on the example above::
     218   
     219        sage: fm = FanMorphism(phi, F1, F2, subdivide=True,
     220        ...                    check=False, verbose=True)
     221        Placing ray images...
     222        Computing chambers...
     223        Subdividing cone 1 of 1...           
     224        sage: fm.domain_fan().is_equivalent(F2)
     225        True
     226    """
     227   
     228    def __init__(self, morphism, domain_fan,
     229                 codomain=None,
     230                 subdivide=False,
     231                 check=True,
     232                 verbose=False):
     233        r"""
     234        Create a fan morphism.
     235       
     236        See :class:`FanMorphism` for documentation.
     237       
     238        TESTS::
     239       
     240            sage: quadrant = Cone([(1,0), (0,1)])
     241            sage: quadrant = Fan([quadrant])
     242            sage: quadrant_bl = quadrant.subdivide([(1,1)])
     243            sage: fm = FanMorphism(identity_matrix(2), quadrant_bl, quadrant)
     244            sage: fm
     245            Fan morphism defined by the matrix
     246            [1 0]
     247            [0 1]
     248            Domain fan: Rational polyhedral fan in 2-d lattice N
     249            Codomain fan: Rational polyhedral fan in 2-d lattice N           
     250            sage: TestSuite(fm).run(skip="_test_category")
     251        """
     252        assert is_Fan(domain_fan)
     253        if is_Fan(codomain):
     254            codomain, codomain_fan = codomain.lattice(), codomain
     255        else:
     256            codomain_fan = None           
     257        if is_FreeModuleMorphism(morphism):
     258            parent = morphism.parent()
     259            A = morphism.matrix()
     260        elif is_Matrix(morphism):
     261            A = morphism
     262            if codomain is None:
     263                raise ValueError("codomain (fan) must be given explicitly if "
     264                                 "morphism is given by a matrix!")
     265            parent = Hom(domain_fan.lattice(), codomain)
     266        else:
     267            raise TypeError("morphism must be either a FreeModuleMorphism "
     268                            "or a matrix!\nGot: %s" % morphism)
     269        super(FanMorphism, self).__init__(parent, A)
     270        self._domain_fan = domain_fan
     271        self._image_cone = dict()
     272        self._preimage_cones = dict()
     273        if codomain_fan is None:
     274            self._construct_codomain_fan()
     275        else:
     276            self._codomain_fan = codomain_fan
     277            if subdivide:
     278                self._subdivide_domain_fan(check, verbose)
     279            elif check:
     280                self._validate()
     281               
     282    def _RISGIS(self):
     283        r"""
     284        Return Ray Images Star Generator Indices Sets.
     285       
     286        OUTPUT:
     287       
     288        - a :class:`tuple` of :class:`frozenset`s of integers, the `i`-th set
     289          is the set of indices of star generators for the minimal cone of the
     290          :meth:`codomain_fan` containing the image of the `i`-th ray of the
     291          :meth:`domain_fan`.
     292         
     293        TESTS::
     294       
     295            sage: diamond = lattice_polytope.octahedron(2)
     296            sage: face = FaceFan(diamond)
     297            sage: normal = NormalFan(diamond)
     298            sage: N = face.lattice()
     299            sage: fm = FanMorphism(identity_matrix(2),
     300            ...           normal, face, subdivide=True)
     301            sage: fm._RISGIS()
     302            (frozenset([3]), frozenset([2]),
     303             frozenset([1]), frozenset([0]),
     304             frozenset([1, 3]), frozenset([0, 1]),
     305             frozenset([0, 2]), frozenset([2, 3]))
     306        """
     307        if "_RISGIS_" not in self.__dict__:
     308            try:
     309                cones = [self._codomain_fan.cone_containing(self(ray))
     310                         for ray in self._domain_fan.rays()]
     311            except ValueError:
     312                self._support_error()
     313            self._RISGIS_ = tuple(frozenset(cone.star_generator_indices())
     314                                  for cone in cones)
     315        return self._RISGIS_
     316           
     317    def _chambers(self):
     318        r"""
     319        Return chambers in the domain corresponding to the codomain fan.
     320       
     321        This function is used during automatic refinement of the domain fans,
     322        see :meth:`_subdivide_domain_fan`.
     323       
     324        OUTPUT:
     325       
     326        - a :class:`tuple` ``(chambers, cone_to_chamber)``, where
     327       
     328        - ``chambers`` is a :class:`list` of :class:`cones
     329          <sage.geometry.cone.ConvexRationalPolyhedralCone>` in the domain of
     330          ``self``;
     331       
     332        - ``cone_to_chamber`` is a :class:`list` of integers, if its `i`-th
     333          element is `j`, then the `j`-th element of ``chambers`` is the
     334          inverse image of the `i`-th generating cone of the codomain fan.
     335                   
     336        TESTS::
     337       
     338            sage: F = NormalFan(lattice_polytope.octahedron(2))
     339            sage: N = F.lattice()
     340            sage: H = End(N)
     341            sage: phi = H([N.0, 0])
     342            sage: fm = FanMorphism(phi, F, F, subdivide=True)
     343            sage: fm._chambers()
     344            ([2-d cone in 2-d lattice N,
     345              1-d cone in 2-d lattice N,
     346              2-d cone in 2-d lattice N],
     347             [0, 1, 2, 1])           
     348        """
     349        kernel_rays = []
     350        for ray in self.kernel().basis():
     351            kernel_rays.append(ray)
     352            kernel_rays.append(-ray)
     353        image_rays = []
     354        for ray in self.image().basis():
     355            image_rays.append(ray)
     356            image_rays.append(-ray)
     357        image = Cone(image_rays)
     358        chambers = []
     359        cone_to_chamber = []
     360        for cone in self._codomain_fan:
     361            chamber = Cone([self.lift(ray) for ray in cone.intersection(image)]
     362                           + kernel_rays, lattice=self.domain())
     363            cone_to_chamber.append(len(chambers))
     364            for i, old_chamber in enumerate(chambers):
     365                if old_chamber.is_equivalent(chamber):
     366                    cone_to_chamber[-1] = i
     367                    break
     368            if cone_to_chamber[-1] == len(chambers):
     369                chambers.append(chamber)
     370        return (chambers, cone_to_chamber)
     371       
     372    def _construct_codomain_fan(self):
     373        r"""
     374        Construct the codomain fan as the image of the domain one.
     375       
     376        .. WARNING::
     377       
     378            This method should be called only during initialization.
     379           
     380        OUTPUT:
     381       
     382        - none, but the codomain fan of ``self`` is set to the constucted fan.
     383       
     384        TESTS::
     385       
     386            sage: quadrant = Cone([(1,0), (0,1)])
     387            sage: quadrant = Fan([quadrant])
     388            sage: quadrant_bl = quadrant.subdivide([(1,1)])
     389            sage: fm = FanMorphism(identity_matrix(2), quadrant_bl)
     390            Traceback (most recent call last):
     391            ...
     392            ValueError: codomain (fan) must be given explicitly
     393            if morphism is given by a matrix!
     394            sage: fm = FanMorphism(identity_matrix(2), quadrant_bl, ZZ^2)
     395            sage: fm.codomain_fan().ray_matrix()  # indirect doctest
     396            [1 0 1]
     397            [0 1 1]
     398        """
     399        # We literally try to construct the image fan and hope that it works.
     400        # If it does not, the fan constructor will raise an exception.
     401        domain_fan = self._domain_fan
     402        self._codomain_fan = Fan(cones=(domain_cone.ambient_ray_indices()
     403                                        for domain_cone in domain_fan),
     404                                 rays=(self(ray) for ray in domain_fan.rays()))
     405   
     406    def _latex_(self):
     407        r"""
     408        Return the `\LaTeX` representation of ``self``.
     409       
     410        OUTPUT:
     411       
     412        - a :class:`string`.
     413       
     414        EXAMPLES::
     415       
     416            sage: quadrant = Cone([(1,0), (0,1)])
     417            sage: quadrant = Fan([quadrant])
     418            sage: quadrant_bl = quadrant.subdivide([(1,1)])
     419            sage: fm = FanMorphism(identity_matrix(2), quadrant_bl, quadrant)
     420            sage: fm
     421            Fan morphism defined by the matrix
     422            [1 0]
     423            [0 1]
     424            Domain fan: Rational polyhedral fan in 2-d lattice N
     425            Codomain fan: Rational polyhedral fan in 2-d lattice N           
     426            sage: print fm._latex_()
     427            \left(\begin{array}{rr}
     428            1 & 0 \\
     429            0 & 1
     430            \end{array}\right) : \Sigma^{2} \to \Sigma^{2}
     431        """
     432        return (r"%s : %s \to %s" % (latex(self.matrix()),
     433                        latex(self.domain_fan()), latex(self.codomain_fan())))
     434           
     435    def _repr_(self):
     436        r"""
     437        Return the string representation of ``self``.
     438       
     439        OUTPUT:
     440       
     441        - a :class:`string`.
     442       
     443        EXAMPLES::
     444       
     445            sage: quadrant = Cone([(1,0), (0,1)])
     446            sage: quadrant = Fan([quadrant])
     447            sage: quadrant_bl = quadrant.subdivide([(1,1)])
     448            sage: fm = FanMorphism(identity_matrix(2), quadrant_bl, quadrant)
     449            sage: print fm._repr_()
     450            Fan morphism defined by the matrix
     451            [1 0]
     452            [0 1]
     453            Domain fan: Rational polyhedral fan in 2-d lattice N
     454            Codomain fan: Rational polyhedral fan in 2-d lattice N           
     455        """
     456        return ("Fan morphism defined by the matrix\n"
     457                "%s\n"
     458                "Domain fan: %s\n"
     459                "Codomain fan: %s"
     460                % (self.matrix(), self.domain_fan(), self.codomain_fan()))
     461   
     462    def _subdivide_domain_fan(self, check, verbose):
     463        r"""
     464        Subdivide the domain fan to make it compatible with the codomain fan.
     465       
     466        .. WARNING::
     467       
     468            This method should be called only during initialization.
     469                   
     470        INPUT:
     471       
     472        - ``check`` -- (default: ``True``) if ``False``, some of the
     473          consistency checks will be omitted, which saves time but can
     474          potentially lead to wrong results. Currently, with
     475          ``check=False`` option there will be no check that ``domain_fan``
     476          maps to ``codomain_fan``;
     477         
     478        - ``verbose`` -- (default: ``False``) if ``True``, some timing
     479          information will be printed in the process.
     480         
     481        OUTPUT:
     482       
     483        - none, but the domain fan of self is replaced with its minimal
     484          refinement, if possible. Otherwise a ``ValueError`` exception is
     485          raised.
     486         
     487        TESTS::
     488       
     489            sage: quadrant = Cone([(1,0), (0,1)])
     490            sage: quadrant = Fan([quadrant])
     491            sage: quadrant_bl = quadrant.subdivide([(1,1)])
     492            sage: fm = FanMorphism(identity_matrix(2), quadrant, quadrant_bl)
     493            Traceback (most recent call last):
     494            ...
     495            ValueError: the image of generating cone #0 of the domain fan
     496            is not contained in a single cone of the codomain fan!
     497            sage: fm = FanMorphism(identity_matrix(2), quadrant,
     498            ...                    quadrant_bl, subdivide=True)
     499            sage: fm.domain_fan().ray_matrix()  # indirect doctest
     500            [1 0 1]
     501            [0 1 1]
     502
     503        Now we demonstrate a more subtle example. We take the first quadrant
     504        as our ``domain_fan``. Then we divide the first quadrant into three
     505        cones, throw away the middle one and take the other two as our
     506        ``codomain_fan``. These fans are incompatible with the identity
     507        lattice morphism since the image of ``domain_fan`` is out of the
     508        support of ``codomain_fan``::
     509       
     510            sage: N = ToricLattice(2)
     511            sage: phi = End(N).identity()
     512            sage: F1 = Fan(cones=[(0,1)], rays=[(1,0), (0,1)])
     513            sage: F2 = Fan(cones=[(0,1), (2,3)],
     514            ...            rays=[(1,0), (2,1), (1,2), (0,1)])
     515            sage: FanMorphism(phi, F1, F2)
     516            Traceback (most recent call last):
     517            ...
     518            ValueError: the image of generating cone #0 of the domain fan
     519            is not contained in a single cone of the codomain fan!
     520            sage: FanMorphism(phi, F1, F2, subdivide=True)
     521            Traceback (most recent call last):
     522            ...
     523            ValueError: morphism defined by
     524            [1 0]
     525            [0 1]
     526            does not map
     527            Rational polyhedral fan in 2-d lattice N
     528            into the support of
     529            Rational polyhedral fan in 2-d lattice N!
     530        """
     531        domain_fan = self._domain_fan
     532        codomain_fan = self._codomain_fan
     533        if verbose:
     534            start = walltime()
     535            print "Placing ray images...",
     536        # Figure out where 1-dimensional cones (i.e. rays) are mapped.
     537        RISGIS = self._RISGIS()
     538        if verbose:
     539            print "%.3f ms" % walltime(start)
     540        # Subdivide cones that require it.
     541        chambers = None # preimages of codomain cones, computed if necessary
     542        new_cones = []
     543        for cone_index, domain_cone in enumerate(domain_fan):
     544            if reduce(operator.and_, (RISGIS[i]
     545                                for i in domain_cone.ambient_ray_indices())):
     546                new_cones.append(domain_cone)
     547                continue
     548            dim = domain_cone.dim()
     549            if chambers is None:
     550                if verbose:
     551                    start = walltime()
     552                    print "Computing chambers...",
     553                chambers, cone_to_chamber = self._chambers()
     554                if verbose:
     555                    print "%.3f ms" % walltime(start)
     556            # Subdivide domain_cone.
     557            if verbose:
     558                start = walltime()
     559                print ("Subdividing cone %d of %d..."
     560                       % (cone_index + 1, domain_fan.ngenerating_cones())),
     561            # Only these chambers intersect domain_cone.
     562            containing_chambers = (cone_to_chamber[j]
     563                for j in reduce(operator.or_, (RISGIS[i]
     564                                for i in domain_cone.ambient_ray_indices())))
     565            # We don't care about chambers of small dimension.
     566            containing_chambers = (chambers[i]
     567                                   for i in set(containing_chambers)
     568                                   if chambers[i].dim() >= dim)
     569            parts = (domain_cone.intersection(chamber)
     570                     for chamber in containing_chambers)
     571            # We cannot leave parts as a generator since we use them twice.
     572            parts = [part for part in parts if part.dim() == dim]
     573            if check:
     574                # Check if the subdivision is complete, i.e. there are no
     575                # missing pieces of domain_cone. To do this, we construct a
     576                # fan from the obtained parts and check that interior points
     577                # of boundary cones of this fan are in the interior of the
     578                # original cone. In any case we know that we are constructing
     579                # a valid fan, so passing check=False to Fan(...) is OK.
     580                if verbose:
     581                    print "%.3f ms" % walltime(start)
     582                    start = walltime()
     583                    print "Checking for missing pieces... ",
     584                cone_subdivision = Fan(parts, check=False)
     585                for cone in cone_subdivision(dim - 1):
     586                    if len(cone.star_generators()) == 1:
     587                        if domain_cone.relative_interior_contains(
     588                                                            sum(cone.rays())):
     589                            self._support_error()
     590            new_cones.extend(parts)
     591            if verbose:
     592                print "%.3f ms" % walltime(start)
     593        if len(new_cones) > domain_fan.ngenerating_cones():
     594            # Construct a new fan keeping old rays in the same order
     595            new_rays = list(domain_fan.rays())
     596            for cone in new_cones:
     597                for ray in cone:
     598                    if ray not in new_rays:
     599                        new_rays.append(ray)
     600            # Replace domain_fan, this is OK since this method must be called
     601            # only during initialization of the FanMorphism.
     602            self._domain_fan = Fan(new_cones, new_rays, domain_fan.lattice(),
     603                                   check=False)
     604            # Also remove RISGIS for the old fan
     605            del self._RISGIS_
     606
     607    def _support_error(self):
     608        r"""
     609        Raise a ``ValueError`` exception due to support incompatibility.
     610       
     611        OUTPUT:
     612       
     613        - none, a ``ValueError`` exception is raised.
     614       
     615        TESTS:
     616       
     617        We deliberately construct an invalid morphism for this demonstration::
     618       
     619            sage: quadrant = Cone([(1,0), (0,1)])
     620            sage: quadrant = Fan([quadrant])
     621            sage: quadrant_bl = quadrant.subdivide([(1,1)])
     622            sage: fm = FanMorphism(identity_matrix(2),
     623            ...             quadrant, quadrant_bl, check=False)
     624           
     625        Now we report that the morphism is invalid::
     626       
     627            sage: fm._support_error()
     628            Traceback (most recent call last):
     629            ...
     630            ValueError: morphism defined by
     631            [1 0]
     632            [0 1]
     633            does not map
     634            Rational polyhedral fan in 2-d lattice N
     635            into the support of
     636            Rational polyhedral fan in 2-d lattice N!           
     637        """
     638        raise ValueError("morphism defined by\n"
     639                         "%s\n"
     640                         "does not map\n"
     641                         "%s\n"
     642                         "into the support of\n"
     643                         "%s!"
     644                    % (self.matrix(), self.domain_fan(), self.codomain_fan()))
     645
     646    def _validate(self):
     647        r"""
     648        Check if ``self`` is indeed compatible with domain and codomain fans.
     649       
     650        OUTPUT:
     651       
     652        - none, but a ``ValueError`` exception is raised if there is a cone of
     653          the domain fan of ``self`` which is not completely contained in a
     654          single cone of the codomain fan of ``self``, or if one of these fans
     655          does not sit in the appropriate lattice.
     656         
     657        EXAMPLES::
     658       
     659            sage: N3 = ToricLattice(3, "N3")
     660            sage: N2 = ToricLattice(2, "N2")
     661            sage: H = Hom(N3, N2)
     662            sage: phi = H([N2.0, N2.1, N2.0])
     663            sage: F1 = Fan(cones=[(0,1,2), (1,2,3)],
     664            ...            rays=[(1,1,1), (1,1,-1), (1,-1,1), (1,-1,-1)],
     665            ...            lattice=N3)
     666            sage: F1.ray_matrix()
     667            [ 1  1  1  1]
     668            [ 1  1 -1 -1]
     669            [ 1 -1  1 -1]
     670            sage: [phi(ray) for ray in F1.rays()]
     671            [N2(2, 1), N2(0, 1), N2(2, -1), N2(0, -1)]
     672            sage: F2 = Fan(cones=[(0,1,2), (1,2,3)],
     673            ...            rays=[(1,1,1), (1,1,-1), (1,2,1), (1,2,-1)],
     674            ...            lattice=N3)
     675            sage: F2.ray_matrix()
     676            [ 1  1  1  1]
     677            [ 1  1  2  2]
     678            [ 1 -1  1 -1]
     679            sage: [phi(ray) for ray in F2.rays()]
     680            [N2(2, 1), N2(0, 1), N2(2, 2), N2(0, 2)]
     681            sage: F3 = Fan(cones=[(0,1), (1,2)],
     682            ...            rays=[(1,0), (2,1), (0,1)],
     683            ...            lattice=N2)
     684            sage: FanMorphism(phi, F2, F3)
     685            Fan morphism defined by the matrix
     686            [1 0]
     687            [0 1]
     688            [1 0]
     689            Domain fan: Rational polyhedral fan in 3-d lattice N3
     690            Codomain fan: Rational polyhedral fan in 2-d lattice N2
     691            sage: FanMorphism(phi, F1, F3)  # indirect doctest
     692            Traceback (most recent call last):
     693            ...
     694            ValueError: morphism defined by
     695            [1 0]
     696            [0 1]
     697            [1 0]
     698            does not map
     699            Rational polyhedral fan in 3-d lattice N3
     700            into the support of
     701            Rational polyhedral fan in 2-d lattice N2!
     702        """
     703        domain_fan = self._domain_fan
     704        if domain_fan.lattice() is not self.domain():
     705            raise ValueError("%s does not sit in %s!"
     706                             % (domain_fan, self.domain()))
     707        codomain_fan = self._codomain_fan
     708        if codomain_fan.lattice() is not self.codomain():
     709            raise ValueError("%s does not sit in %s!"
     710                             % (codomain_fan, self.codomain()))
     711        RISGIS = self._RISGIS()
     712        for n, domain_cone in enumerate(domain_fan):
     713            if not reduce(operator.and_,
     714                (RISGIS[i] for i in domain_cone.ambient_ray_indices())):
     715                raise ValueError("the image of generating cone #%d of the "
     716                                 "domain fan is not contained in a single "
     717                                 "cone of the codomain fan!" % n)
     718
     719    def codomain_fan(self):
     720        r"""
     721        Return the codomain fan of ``self``.
     722       
     723        OUTPUT:
     724       
     725        - a :class:`fan <sage.geometry.fan.RationalPolyhedralFan>`.
     726       
     727        EXAMPLES::
     728       
     729            sage: quadrant = Cone([(1,0), (0,1)])
     730            sage: quadrant = Fan([quadrant])
     731            sage: quadrant_bl = quadrant.subdivide([(1,1)])
     732            sage: fm = FanMorphism(identity_matrix(2), quadrant_bl, quadrant)
     733            sage: fm.codomain_fan()
     734            Rational polyhedral fan in 2-d lattice N
     735            sage: fm.codomain_fan() is quadrant
     736            True
     737        """
     738        return self._codomain_fan
     739
     740    def domain_fan(self):
     741        r"""
     742        Return the codomain fan of ``self``.
     743       
     744        OUTPUT:
     745       
     746        - a :class:`fan <sage.geometry.fan.RationalPolyhedralFan>`.
     747       
     748        EXAMPLES::
     749       
     750            sage: quadrant = Cone([(1,0), (0,1)])
     751            sage: quadrant = Fan([quadrant])
     752            sage: quadrant_bl = quadrant.subdivide([(1,1)])
     753            sage: fm = FanMorphism(identity_matrix(2), quadrant_bl, quadrant)
     754            sage: fm.domain_fan()
     755            Rational polyhedral fan in 2-d lattice N
     756            sage: fm.domain_fan() is quadrant_bl
     757            True
     758        """
     759        return self._domain_fan
     760       
     761    def image_cone(self, cone):
     762        r"""
     763        Return the cone of the codomain fan containing the image of ``cone``.
     764       
     765        INPUT:
     766       
     767        - ``cone`` -- a :class:`cone
     768          <sage.geometry.cone.ConvexRationalPolyhedralCone>` equivalent to a
     769          cone of the :meth:`domain_fan` of ``self``.
     770         
     771        OUTPUT:
     772       
     773        - a :class:`cone <sage.geometry.cone.ConvexRationalPolyhedralCone>`
     774          of the :meth:`codomain_fan` of ``self``.
     775         
     776        EXAMPLES::
     777       
     778            sage: quadrant = Cone([(1,0), (0,1)])
     779            sage: quadrant = Fan([quadrant])
     780            sage: quadrant_bl = quadrant.subdivide([(1,1)])
     781            sage: fm = FanMorphism(identity_matrix(2), quadrant_bl, quadrant)
     782            sage: fm.image_cone(Cone([(1,0)]))
     783            1-d cone of Rational polyhedral fan in 2-d lattice N
     784            sage: fm.image_cone(Cone([(1,1)]))
     785            2-d cone of Rational polyhedral fan in 2-d lattice N
     786           
     787        TESTS:
     788       
     789        We check that complete codomain fans are handled correctly, since a
     790        different algorithm is used in this case::
     791       
     792            sage: diamond = lattice_polytope.octahedron(2)
     793            sage: face = FaceFan(diamond)
     794            sage: normal = NormalFan(diamond)
     795            sage: N = face.lattice()
     796            sage: fm = FanMorphism(identity_matrix(2),
     797            ...           normal, face, subdivide=True)
     798            sage: fm.image_cone(Cone([(1,0)]))
     799            1-d cone of Rational polyhedral fan in 2-d lattice N
     800            sage: fm.image_cone(Cone([(1,1)]))
     801            2-d cone of Rational polyhedral fan in 2-d lattice N
     802        """
     803        cone = self._domain_fan.embed(cone)
     804        if cone not in self._image_cone:
     805            codomain_fan = self._codomain_fan()
     806            if cone.is_trivial():
     807                self._image_cone[cone] = codomain_fan(0)[0]
     808            elif codomain_fan.is_complete():
     809                # Optimization for a common case
     810                RISGIS = self._RISGIS()
     811                CSGIS = set(reduce(operator.and_,
     812                            (RISGIS[i] for i in cone.ambient_ray_indices())))
     813                image_cone = codomain_fan.generating_cone(CSGIS.pop())
     814                for i in CSGIS:
     815                    image_cone = image_cone.intersection(
     816                                            codomain_fan.generating_cone(i))
     817                self._image_cone[cone] = image_cone
     818            else:
     819                self._image_cone[cone] = codomain_fan.cone_containing(
     820                                                    self(ray) for ray in cone)
     821        return self._image_cone[cone]
     822
     823    def preimage_cones(self, cone):
     824        r"""
     825        Return cones of the domain fan whose :meth:`image_cone` is ``cone``.
     826       
     827        INPUT:
     828       
     829        - ``cone`` -- a :class:`cone
     830          <sage.geometry.cone.ConvexRationalPolyhedralCone>` equivalent to a
     831          cone of the :meth:`codomain_fan` of ``self``.
     832         
     833        OUTPUT:
     834       
     835        - a :class:`tuple` of :class:`cones
     836          <sage.geometry.cone.ConvexRationalPolyhedralCone>` of the
     837          :meth:`domain_fan` of ``self``.
     838         
     839        EXAMPLES::
     840       
     841            sage: quadrant = Cone([(1,0), (0,1)])
     842            sage: quadrant = Fan([quadrant])
     843            sage: quadrant_bl = quadrant.subdivide([(1,1)])
     844            sage: fm = FanMorphism(identity_matrix(2), quadrant_bl, quadrant)
     845            sage: fm.preimage_cones(Cone([(1,0)]))
     846            (1-d cone of Rational polyhedral fan in 2-d lattice N,)
     847            sage: fm.preimage_cones(Cone([(1,0), (0,1)]))
     848            (1-d cone of Rational polyhedral fan in 2-d lattice N,
     849             2-d cone of Rational polyhedral fan in 2-d lattice N,
     850             2-d cone of Rational polyhedral fan in 2-d lattice N)
     851             
     852        TESTS:
     853           
     854        We check that reviewer's example from Trac #9972 is handled correctly::
     855       
     856            sage: N1 = ToricLattice(1)
     857            sage: N2 = ToricLattice(2)
     858            sage: Hom21 = Hom(N2, N1)
     859            sage: pr = Hom21([N1.0,0])
     860            sage: P1xP1 = toric_varieties.P1xP1()
     861            sage: f = FanMorphism(pr, P1xP1.fan())
     862            sage: c = f.image_cone(Cone([(1,0), (0,1)]))
     863            sage: c
     864            1-d cone of Rational polyhedral fan in 1-d lattice N
     865            sage: f.preimage_cones(c)
     866            (1-d cone of Rational polyhedral fan in 2-d lattice N,
     867             2-d cone of Rational polyhedral fan in 2-d lattice N,
     868             2-d cone of Rational polyhedral fan in 2-d lattice N)
     869        """
     870        cone = self._codomain_fan.embed(cone)
     871        domain_fan = self._domain_fan
     872        if cone not in self._preimage_cones:
     873            # "Images of preimages" must sit in all of these generating cones
     874            CSGI = cone.star_generator_indices()
     875            RISGIS = self._RISGIS()
     876            possible_rays = frozenset(i for i in range(domain_fan.nrays())
     877                                        if RISGIS[i].issuperset(CSGI))
     878            preimage_cones = []
     879            for dcones in domain_fan.cones():
     880                for dcone in dcones:
     881                    if (possible_rays.issuperset(dcone.ambient_ray_indices())
     882                        and self.image_cone(dcone) == cone):
     883                        preimage_cones.append(dcone)
     884            self._preimage_cones[cone] = tuple(preimage_cones)
     885        return self._preimage_cones[cone]