Ticket #13517: trac_#13517_new_class_reviewed.patch

File trac_#13517_new_class_reviewed.patch, 8.1 KB (added by moritz, 9 years ago)

reviewed

  • sage/geometry/all.py

    # HG changeset patch
    # User Moritz Firsching <moritz@math.fu-berlin.de>
    # Date 1353607317 -3600
    # Node ID bbda54ee32f78194100428f35f7646670365b25d
    # Parent  29442626e0913ea9f52957230459e89fa10708a4
    trac #13517: reviewed version with arbitrary base_ring
    
    diff --git a/sage/geometry/all.py b/sage/geometry/all.py
    a b  
    1717
    1818
    1919import toric_plotter
     20
     21from voronoi_diagram import *
  • sage/geometry/all.py

    # HG changeset patch
    # User Moritz Firsching <moritz@math.fu-berlin.de>
    # Date 1353607937 -3600
    # Node ID 0c621a37a1a6c9a3b351d537166fb5d89661e50a
    # Parent  29442626e0913ea9f52957230459e89fa10708a4
    Trac: #13517 reviewed first version with arbitrary base_ring
    
    diff --git a/sage/geometry/all.py b/sage/geometry/all.py
    a b  
    1717
    1818
    1919import toric_plotter
     20
     21from voronoi_diagram import *
  • new file sage/geometry/voronoi_diagram.py

    diff --git a/sage/geometry/voronoi_diagram.py b/sage/geometry/voronoi_diagram.py
    new file mode 100644
    - +  
     1r"""
     2Voronoi diagram
     3
     4This module provides the class :class: `voronoi_diagram` for computing the
     5Voronoi diagram of a finite list of points in \RR^d.
     6"""
     7
     8#*****************************************************************************
     9#       Copyright (C) 2012 Moritz Firsching <moritz@math.fu-berlin.de>
     10#
     11#  Distributed under the terms of the GNU General Public License (GPL)
     12#
     13#                  http://www.gnu.org/licenses/
     14#*****************************************************************************
     15
     16from sage.structure.sage_object import SageObject
     17from sage.geometry.polyhedron.constructor import Polyhedron
     18from sage.all import RDF, QQ
     19from sage.geometry.triangulation.point_configuration import PointConfiguration
     20from sage.modules.all import vector
     21from sage.plot.all import line, point, rainbow, plot
     22
     23
     24class voronoi_diagram(SageObject):
     25    r"""
     26    base class for the  Voronoi diagram.
     27    Computes the Voronoi diagram of a list of points
     28
     29    INPUT:
     30
     31    - ``points`` a list of points. Any valid input for the
     32    :class:`PointConfiguration` will do.
     33
     34    EXAMPLES:
     35        Get the Voronoi diagram for some points in \RR^3
     36        ::
     37        sage: V=voronoi_diagram([[1,3,.3],[2,-2,1],[-1,2,-.1]]); V
     38        The Voronoi diagram of 3 points.           
     39        sage: voronoi_diagram([])
     40        The empty Voronoi diagram.
     41       
     42    ALGORITHM:
     43
     44    We use hyperplanes tangent to the paraboloid one dimension higher to
     45    get a convex polyhedron and then project back to one dimension lower
     46
     47    See for example [M2002]
     48
     49    REFERENCES:
     50
     51    ..  [M2002]
     52        Jiri Matousek,
     53        "Lectures on Discrete Geometry", Springer, Ch.5.7, p.118.
     54
     55    AUTHORS:
     56
     57    - Moritz Firsching (2012-09-21)
     58    """
     59    def __init__(self, points):
     60        r"""
     61        See ``voronoi_diagram`` for full documentation.
     62        EXAMPLES::
     63        sage: V=voronoi_diagram([[1,3,3],[2,-2,1],[-1,2,-1]]); V
     64        The Voronoi diagram of 3 points. 
     65       
     66        """
     67        self._P=[]
     68        self._points=PointConfiguration(points)
     69        self._n=self._points.n_points()
     70        if self._n==0 or self._points.base_ring().is_subring(QQ):
     71            self._base_ring=QQ
     72        else:
     73            self._base_ring=RDF
     74           
     75        if self._n:
     76            self._d=self._points.ambient_dim()
     77            e=[([sum(vector(i)[k]**2 for k in
     78            range(self._d))]+[(-2)*vector(i)[l] for l in range(self._d)]+[1])
     79            for i in self._points]
     80            e=[[self._base_ring(i) for i in k] for k in e]
     81            p=Polyhedron(ieqs = e, base_ring=self._base_ring)
     82        for i in range(self._n):
     83            equ=p.Hrepresentation(i)
     84            pvert=[[u[k] for k in range(self._d)] for u in equ.incident() if
     85            u.is_vertex()]
     86            prays=[[u[k] for k in range(self._d)] for u in equ.incident() if
     87            u.is_ray()]
     88            pline=[[u[k] for k in range(self._d)] for u in equ.incident() if
     89            u.is_line()]
     90            (self._P).append(Polyhedron(vertices=pvert, lines=pline, rays=prays,
     91             base_ring=self._base_ring))
     92
     93    def points(self):
     94        r""" Returns the input points (as a PointConfiguration).
     95
     96        EXAMPLES::
     97
     98            sage: V=voronoi_diagram([[.5,3],[2,5],[4,5],[4,-1]]); V.points()
     99            A point configuration in QQ^2 consisting of 4 points. The
     100            triangulations of this point configuration are assumed to be
     101            connected, not necessarily
     102            fine, not necessarily regular.
     103        """
     104        return self._points
     105
     106    def ambient_dim(self):
     107        r"""
     108        Returns the ambient dimension of the points.
     109
     110        EXAMPLES::
     111            sage: V=voronoi_diagram([[.5,3],[2,5],[4,5],[4,-1]])
     112            sage: V.ambient_dim()
     113            2
     114            sage: V=voronoi_diagram([[1,2,3,4,5,6]]); V.ambient_dim()
     115            6
     116        """
     117        return self._d
     118
     119    def regions(self):
     120        r"""
     121        Returns the Voronoi regions of the Voronoi diagram as a list of
     122        polyhedra.
     123
     124        EXAMPLES::
     125            sage: V=voronoi_diagram([[1,3,.3],[2,-2,1],[-1,2,-.1]]); V.regions()
     126            [A 3-dimensional polyhedron in RDF^3 defined as the convex hull of 1
     127             vertex, 2 rays, 1 line, A 3-dimensional polyhedron in RDF^3 defined
     128              as the
     129            convex hull of 1 vertex, 2 rays, 1 line, A 3-dimensional polyhedron
     130            in RDF^3 defined
     131            as the convex hull of 1 vertex, 2 rays, 1 line]
     132        """
     133        return self._P
     134
     135    def base_ring(self):
     136        r"""
     137        Returns the base_ring of the regions of the Voronoi diagram.
     138
     139        EXAMPLES::
     140            sage: V=voronoi_diagram([[1,3,1],[2,-2,1],[-1,2,1/2]]); V.base_ring()
     141            Rational Field
     142            sage: V=voronoi_diagram([[1,3.14],[2,-2/3],[-1,22]]); V.base_ring()
     143            Real Double Field
     144            sage: V=voronoi_diagram([[1,3],[2,4]]); V.base_ring()
     145            Rational Field
     146        """
     147        return self._base_ring
     148
     149    def _repr_(self):
     150        r"""
     151        Return a description of the Voronoi diagram.
     152
     153        EXAMPLES::
     154            sage: V=voronoi_diagram([[1,3,.3],[2,-2,1],[-1,2,-.1]]); V
     155            The Voronoi diagram of 3 points.
     156            sage: voronoi_diagram([])
     157            The empty Voronoi diagram.
     158
     159        """
     160        desc = ''
     161        if self._n:
     162            desc+= 'The Voronoi diagram of '
     163            desc+=str(self._n)
     164            desc+= ' points.'
     165        else:
     166            desc+='The empty Voronoi diagram.'
     167
     168        return desc
     169
     170
     171    def plot(self, **kwds):
     172        """
     173        Return a graphical representation for 2-dimensional Voronoi diagrams.
     174
     175        INPUT:
     176
     177        - ``**kwds`` -- optional keyword parameters, passed on as arguments for
     178        plot().
     179
     180        OUTPUT:
     181
     182        A graphics object.
     183
     184        EXAMPLES::
     185            sage: P=[[0.671, 0.650],[0.258, 0.767], [0.562, 0.406],\
     186            [0.254, 0.709], [0.493, 0.879]]
     187
     188            sage: V=voronoi_diagram(P); S=V.plot()
     189
     190            sage: show(S,xmin=0,xmax=1,ymin=0,ymax=1, aspect_ratio=1,\
     191            axes=false)
     192
     193        Trying to plot a Voronoi diagram of dimension other than 2 gives an
     194        error::
     195            sage: voronoi_diagram([[1,2,3],[6,5,4]]).plot()
     196            Traceback (most recent call last):
     197            ...
     198            NotImplementedError: Plotting of 3-dimensional Voronoi diagrams not
     199                implemented
     200
     201        """
     202        if self.ambient_dim()==2:
     203            S=line([])
     204            for i,j in enumerate(self._points):
     205                S+=(self._P[i]).render_solid(color=rainbow(self._n)[i],
     206                zorder=1)
     207                S+=point(j, color=rainbow(self._n)[i], pointsize=10,zorder=3)
     208                S+=point(vector(j), color='black',pointsize=20,zorder=2)
     209            return plot(S,**kwds)
     210        raise NotImplementedError('Plotting of '+str(self.ambient_dim())+
     211                                  '-dimensional Voronoi diagrams not'+
     212                                  ' implemented')