Ticket #11763: trac_11763_ZZ_polyhedron-rebase.patch
File trac_11763_ZZ_polyhedron-rebase.patch, 14.1 KB (added by , 8 years ago) |
---|
-
doc/en/reference/geometry.rst
# HG changeset patch # User Volker Braun <vbraun@stp.dias.ie> # Date 1315501129 14400 # Node ID 01ee8583d6d33ad70ced4ac52f0c5c237335fddb # Parent fcf92c73da9d7fd61a4cdc525e877a7c9a9f603e Trac #11763: Add ZZ as allowed base ring for polyhedra diff --git a/doc/en/reference/geometry.rst b/doc/en/reference/geometry.rst
a b 25 25 sage/geometry/polyhedron/plot 26 26 sage/geometry/polyhedron/base 27 27 sage/geometry/polyhedron/base_QQ 28 sage/geometry/polyhedron/base_ZZ 28 29 sage/geometry/polyhedron/base_RDF 29 30 sage/geometry/polyhedron/backend_cdd 30 31 sage/geometry/polyhedron/backend_ppl -
sage/geometry/polyhedron/backend_ppl.py
diff --git a/sage/geometry/polyhedron/backend_ppl.py b/sage/geometry/polyhedron/backend_ppl.py
a b 11 11 Variable, Linear_Expression, 12 12 line, ray, point ) 13 13 14 from base import Polyhedron_base 14 15 from base_QQ import Polyhedron_QQ 16 from base_ZZ import Polyhedron_ZZ 15 17 from representation import ( 16 18 PolyhedronRepresentation, 17 19 Hrepresentation, … … 22 24 23 25 24 26 ######################################################################### 25 class Polyhedron_ QQ_ppl(Polyhedron_QQ):27 class Polyhedron_ppl(Polyhedron_base): 26 28 """ 27 Polyhedra over `\QQ`with ppl29 Polyhedra with ppl 28 30 29 31 INPUT: 30 32 … … 63 65 EXAMPLES:: 64 66 65 67 sage: p = Polyhedron(backend='ppl') 66 sage: from sage.geometry.polyhedron.backend_ppl import Polyhedron_ QQ_ppl67 sage: Polyhedron_ QQ_ppl._init_from_Vrepresentation(p, 2, [], [], [])68 sage: from sage.geometry.polyhedron.backend_ppl import Polyhedron_ppl 69 sage: Polyhedron_ppl._init_from_Vrepresentation(p, 2, [], [], []) 68 70 """ 69 71 gs = Generator_System() 70 72 if vertices is None: vertices = [] … … 106 108 EXAMPLES:: 107 109 108 110 sage: p = Polyhedron(backend='ppl') 109 sage: from sage.geometry.polyhedron.backend_ppl import Polyhedron_ QQ_ppl110 sage: Polyhedron_ QQ_ppl._init_from_Hrepresentation(p, 2, [], [])111 sage: from sage.geometry.polyhedron.backend_ppl import Polyhedron_ppl 112 sage: Polyhedron_ppl._init_from_Hrepresentation(p, 2, [], []) 111 113 """ 112 114 cs = Constraint_System() 113 115 if ieqs is None: ieqs = [] … … 212 214 The empty polyhedron in QQ^0 213 215 sage: Polyhedron(backend='ppl')._init_empty_polyhedron(0) 214 216 """ 215 super(Polyhedron_ QQ_ppl, self)._init_empty_polyhedron(ambient_dim)217 super(Polyhedron_ppl, self)._init_empty_polyhedron(ambient_dim) 216 218 self._ppl_polyhedron = C_Polyhedron(ambient_dim, 'empty') 217 219 218 220 221 222 223 ######################################################################### 224 class Polyhedron_QQ_ppl(Polyhedron_ppl, Polyhedron_QQ): 225 """ 226 Polyhedra over `\QQ` with ppl 227 228 INPUT: 229 230 - ``ambient_dim`` -- integer. The dimension of the ambient space. 231 232 - ``Vrep`` -- a list ``[vertices, rays, lines]``. 233 234 - ``Hrep`` -- a list ``[ieqs, eqns]``. 235 236 EXAMPLES:: 237 238 sage: p = Polyhedron(vertices=[(0,0),(1,0),(0,1)], rays=[(1,1)], lines=[], 239 ... backend='ppl', base_ring=QQ) 240 sage: TestSuite(p).run(skip='_test_pickling') 241 """ 242 pass 243 244 245 ######################################################################### 246 class Polyhedron_ZZ_ppl(Polyhedron_ppl, Polyhedron_ZZ): 247 """ 248 Polyhedra over `\ZZ` with ppl 249 250 INPUT: 251 252 - ``ambient_dim`` -- integer. The dimension of the ambient space. 253 254 - ``Vrep`` -- a list ``[vertices, rays, lines]``. 255 256 - ``Hrep`` -- a list ``[ieqs, eqns]``. 257 258 EXAMPLES:: 259 260 sage: p = Polyhedron(vertices=[(0,0),(1,0),(0,1)], rays=[(1,1)], lines=[]) 261 ... backend='ppl', base_ring=ZZ) 262 sage: TestSuite(p).run(skip='_test_pickling') 263 """ 264 pass -
sage/geometry/polyhedron/base.py
diff --git a/sage/geometry/polyhedron/base.py b/sage/geometry/polyhedron/base.py
a b 537 537 desc += 'A ' + repr(self.dim()) + '-dimensional polyhedron' 538 538 desc += ' in ' 539 539 if self.field()==QQ: desc += 'QQ' 540 else: desc += 'RDF' 540 if self.base_ring() is QQ: desc += 'QQ' 541 elif self.base_ring() is ZZ: desc += 'ZZ' 542 elif self.base_ring() is RDF: desc += 'RDF' 543 else: assert False 541 544 desc += '^' + repr(self.ambient_dim()) 542 545 543 546 if self.n_vertices()>0: … … 1350 1353 sage: poly_test.ambient_space() 1351 1354 Vector space of dimension 4 over Rational Field 1352 1355 """ 1353 from sage.modules.free_module import VectorSpace1354 return VectorSpace(self.base_ring(), self.ambient_dim())1356 from sage.modules.free_module import FreeModule 1357 return FreeModule(self.base_ring(), self.ambient_dim()) 1355 1358 1356 1359 Vrepresentation_space = ambient_space 1357 1360 -
new file sage/geometry/polyhedron/base_ZZ.py
diff --git a/sage/geometry/polyhedron/base_ZZ.py b/sage/geometry/polyhedron/base_ZZ.py new file mode 100644
- + 1 """ 2 Base class for polyhedra over `\ZZ` 3 """ 4 5 ######################################################################## 6 # Copyright (C) 2011 Volker Braun <vbraun.name@gmail.com> 7 # 8 # Distributed under the terms of the GNU General Public License (GPL) 9 # 10 # http://www.gnu.org/licenses/ 11 ######################################################################## 12 13 14 15 from sage.rings.all import ZZ, QQ 16 from sage.misc.all import cached_method 17 from sage.matrix.constructor import matrix 18 19 from constructor import Polyhedron 20 from base import Polyhedron_base 21 22 23 24 ######################################################################### 25 class Polyhedron_ZZ(Polyhedron_base): 26 """ 27 Base class for Polyhedra over `\ZZ` 28 29 TESTS:: 30 31 sage: p = Polyhedron([(0,0)], base_ring=ZZ); p 32 A 0-dimensional polyhedron in ZZ^2 defined as the convex hull of 1 vertex 33 sage: TestSuite(p).run(skip='_test_pickling') 34 """ 35 def _is_zero(self, x): 36 """ 37 Test whether ``x`` is zero. 38 39 INPUT: 40 41 - ``x`` -- a number in the base ring. 42 43 OUTPUT: 44 45 Boolean. 46 47 EXAMPLES:: 48 49 sage: p = Polyhedron([(0,0)], base_ring=ZZ) 50 sage: p._is_zero(0) 51 True 52 sage: p._is_zero(1/100000) 53 False 54 """ 55 return x==0 56 57 def _is_nonneg(self, x): 58 """ 59 Test whether ``x`` is nonnegative. 60 61 INPUT: 62 63 - ``x`` -- a number in the base ring. 64 65 OUTPUT: 66 67 Boolean. 68 69 EXAMPLES:: 70 71 sage: p = Polyhedron([(0,0)], base_ring=ZZ) 72 sage: p._is_nonneg(1) 73 True 74 sage: p._is_nonneg(-1/100000) 75 False 76 """ 77 return x>=0 78 79 def _is_positive(self, x): 80 """ 81 Test whether ``x`` is positive. 82 83 INPUT: 84 85 - ``x`` -- a number in the base ring. 86 87 OUTPUT: 88 89 Boolean. 90 91 EXAMPLES:: 92 93 sage: p = Polyhedron([(0,0)], base_ring=ZZ) 94 sage: p._is_positive(1) 95 True 96 sage: p._is_positive(0) 97 False 98 """ 99 return x>0 100 101 _base_ring = ZZ 102 103 def is_lattice_polytope(self): 104 r""" 105 Return whether the polyhedron is a lattice polytope. 106 107 OUTPUT: 108 109 ``True`` if the polyhedron is compact and has only integral 110 vertices, ``False`` otherwise. 111 112 EXAMPLES:: 113 114 sage: polytopes.cross_polytope(3).is_lattice_polytope() 115 True 116 sage: polytopes.regular_polygon(5).is_lattice_polytope() 117 False 118 """ 119 return True 120 121 @cached_method 122 def polar(self): 123 """ 124 Return the polar (dual) polytope. 125 126 The polytope must have the IP-property (see 127 :meth:`has_IP_property`), that is, the origin must be an 128 interior point. In particular, it must be full-dimensional. 129 130 OUTPUT: 131 132 The polytope whose vertices are the coefficient vectors of the 133 inequalities of ``self`` with inhomogeneous term normalized to 134 unity. 135 136 EXAMPLES:: 137 138 sage: p = Polyhedron(vertices=[(1,0,0),(0,1,0),(0,0,1),(-1,-1,-1)], base_ring=ZZ) 139 sage: p.polar() 140 A 3-dimensional polyhedron in ZZ^3 defined as the convex hull of 4 vertices 141 sage: type(_) 142 <class 'sage.geometry.polyhedron.backend_ppl.Polyhedron_ZZ_ppl'> 143 sage: p.polar().base_ring() 144 Integer Ring 145 """ 146 if not self.has_IP_property(): 147 raise ValueError('The polytope must have the IP property.') 148 149 vertices = [ ieq.A()/ieq.b() for 150 ieq in self.inequality_generator() ] 151 if all( all(v_i in ZZ for v_i in v) for v in vertices): 152 return Polyhedron(vertices=vertices, base_ring=ZZ) 153 else: 154 return Polyhedron(vertices=vertices, base_ring=QQ) 155 156 @cached_method 157 def is_reflexive(self): 158 """ 159 EXAMPLES:: 160 161 sage: p = Polyhedron(vertices=[(1,0,0),(0,1,0),(0,0,1),(-1,-1,-1)], base_ring=ZZ) 162 sage: p.is_reflexive() 163 True 164 """ 165 return self.polar().is_lattice_polytope() 166 167 @cached_method 168 def has_IP_property(self): 169 """ 170 Test whether the polyhedron has the IP property. 171 172 The IP (interior point) property means that 173 174 * ``self`` is compact (a polytope). 175 176 * ``self`` contains the origin. 177 178 This implies that 179 180 * ``self`` is full-dimensional. 181 182 * The dual polyhedron is again a polytope (that is, a compact 183 polyhedron), though not necessarily a lattice polytope. 184 185 REFERENCES:: 186 187 .. [PALP] 188 Maximilian Kreuzer, Harald Skarke: 189 "PALP: A Package for Analyzing Lattice Polytopes 190 with Applications to Toric Geometry" 191 Comput.Phys.Commun. 157 (2004) 87-106 192 http://arxiv.org/abs/math/0204356 193 """ 194 return self.is_compact() and self.contains(self.ambient_space().zero()) 195 196 def fibrations(self): 197 """ 198 EXAMPLES:: 199 200 sage: P = Polyhedron(toric_varieties.P4_11169().fan().rays(), base_ring=ZZ) 201 sage: list( P.fibrations() ) 202 [A 2-dimensional polyhedron in QQ^4 defined as the convex hull of 3 vertices] 203 """ 204 if not self.has_IP_property(): 205 raise ValueError('The polytope must have the IP property.') 206 # todo: just go through pairs of points on the same facet 207 fiber_dimension = 2 208 points = self.integral_points() 209 fibers = set() 210 for i1 in range(0,len(points)): 211 p1 = points[i1] 212 if p1.is_zero(): 213 continue 214 for i2 in range(i1+1,len(points)): 215 p2 = points[i2] 216 if p2.is_zero(): 217 continue 218 plane_12 = Polyhedron(lines=[p1,p2]) 219 if plane_12.dim() != 2: 220 continue 221 fiber = self.intersection(plane_12) 222 if not fiber.is_lattice_polytope(): 223 continue 224 fiber_matrix = matrix(ZZ,sorted(fiber.integral_points())) 225 fiber_matrix.set_immutable() 226 if fiber_matrix not in fibers: 227 yield fiber 228 fibers.update([fiber_matrix]) 229 230 231 232 233 234 235 236 237 -
sage/geometry/polyhedron/constructor.py
diff --git a/sage/geometry/polyhedron/constructor.py b/sage/geometry/polyhedron/constructor.py
a b 159 159 * ``'cddf'``: cdd with floating-point coefficients 160 160 161 161 * ``'ppl'``: use ppl 162 (:mod:`~sage.geometry.polyhedron.backend_ppl`) with `\QQ` 163 coefficients.162 (:mod:`~sage.geometry.polyhedron.backend_ppl`) with `\QQ` or 163 `\ZZ` coefficients depending on ``base_ring``. 164 164 165 165 Some backends support further optional arguments: 166 166 … … 272 272 ambient_dim = 0 273 273 274 274 if backend is not None: 275 if backend=='ppl' :275 if backend=='ppl' and base_ring is QQ: 276 276 from backend_ppl import Polyhedron_QQ_ppl 277 277 return Polyhedron_QQ_ppl(ambient_dim, Vrep, Hrep, minimize=minimize) 278 if backend=='ppl' and base_ring is ZZ: 279 from backend_ppl import Polyhedron_ZZ_ppl 280 return Polyhedron_ZZ_ppl(ambient_dim, Vrep, Hrep, minimize=minimize) 278 281 if backend=='cddr': 279 282 from backend_cdd import Polyhedron_QQ_cdd 280 283 return Polyhedron_QQ_cdd(ambient_dim, Vrep, Hrep, verbose=verbose) 281 284 if backend=='cddf': 282 285 from backend_cdd import Polyhedron_RDF_cdd 283 286 return Polyhedron_RDF_cdd(ambient_dim, Vrep, Hrep, verbose=verbose) 287 raise ValueError('No suitable backend implemented.') 284 288 285 289 if base_ring is QQ: 286 290 from backend_ppl import Polyhedron_QQ_ppl … … 288 292 elif base_ring is RDF: 289 293 from backend_cdd import Polyhedron_RDF_cdd 290 294 return Polyhedron_RDF_cdd(ambient_dim, Vrep, Hrep, verbose=verbose) 295 elif base_ring is ZZ: 296 from backend_ppl import Polyhedron_ZZ_ppl 297 return Polyhedron_ZZ_ppl(ambient_dim, Vrep, Hrep, minimize=minimize) 291 298 else: 292 299 raise ValueError('Polyhedron objects can only be constructed over QQ and RDF') 293 300 -
sage/geometry/polyhedron/representation.py
diff --git a/sage/geometry/polyhedron/representation.py b/sage/geometry/polyhedron/representation.py
a b 60 60 (1, 2, 3) 61 61 sage: TestSuite(pr).run(skip='_test_pickling') 62 62 """ 63 self._representation_data = tuple(data) 63 R = polyhedron.base_ring() 64 self._representation_data = tuple(R(x) for x in data) 64 65 self._polyhedron = polyhedron 65 66 66 67 def __len__(self):