Fix implementation of Polyhedron.__new__()
[linpy.git] / linpy / domains.py
1 # Copyright 2014 MINES ParisTech
2 #
3 # This file is part of LinPy.
4 #
5 # LinPy is free software: you can redistribute it and/or modify
6 # it under the terms of the GNU General Public License as published by
7 # the Free Software Foundation, either version 3 of the License, or
8 # (at your option) any later version.
9 #
10 # LinPy is distributed in the hope that it will be useful,
11 # but WITHOUT ANY WARRANTY; without even the implied warranty of
12 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 # GNU General Public License for more details.
14 #
15 # You should have received a copy of the GNU General Public License
16 # along with LinPy. If not, see <http://www.gnu.org/licenses/>.
17
18 import ast
19 import functools
20 import re
21 import math
22
23 from fractions import Fraction
24
25 from . import islhelper
26 from .islhelper import mainctx, libisl
27 from .linexprs import LinExpr, Symbol, Rational
28 from .geometry import GeometricObject, Point, Vector
29
30
31 __all__ = [
32 'Domain',
33 'And', 'Or', 'Not',
34 ]
35
36
37 @functools.total_ordering
38 class Domain(GeometricObject):
39 """
40 A domain is a union of polyhedra. Unlike polyhedra, domains allow exact
41 computation of union and complementary operations.
42
43 A domain with a unique polyhedron is automatically subclassed as a
44 Polyhedron instance.
45 """
46
47 __slots__ = (
48 '_polyhedra',
49 '_symbols',
50 '_dimension',
51 )
52
53 def __new__(cls, *polyhedra):
54 """
55 Return a domain from a sequence of polyhedra.
56
57 >>> square = Polyhedron('0 <= x <= 2, 0 <= y <= 2')
58 >>> square2 = Polyhedron('2 <= x <= 4, 2 <= y <= 4')
59 >>> dom = Domain([square, square2])
60
61 It is also possible to build domains from polyhedra using arithmetic
62 operators Domain.__and__(), Domain.__or__() or functions And() and Or(),
63 using one of the following instructions:
64
65 >>> square = Polyhedron('0 <= x <= 2, 0 <= y <= 2')
66 >>> square2 = Polyhedron('2 <= x <= 4, 2 <= y <= 4')
67 >>> dom = square | square2
68 >>> dom = Or(square, square2)
69
70 Alternatively, a domain can be built from a string:
71
72 >>> dom = Domain('0 <= x <= 2, 0 <= y <= 2; 2 <= x <= 4, 2 <= y <= 4')
73
74 Finally, a domain can be built from a GeometricObject instance, calling
75 the GeometricObject.asdomain() method.
76 """
77 from .polyhedra import Polyhedron
78 if len(polyhedra) == 1:
79 argument = polyhedra[0]
80 if isinstance(argument, str):
81 return cls.fromstring(argument)
82 elif isinstance(argument, GeometricObject):
83 return argument.aspolyhedron()
84 else:
85 raise TypeError('argument must be a string '
86 'or a GeometricObject instance')
87 else:
88 for polyhedron in polyhedra:
89 if not isinstance(polyhedron, Polyhedron):
90 raise TypeError('arguments must be Polyhedron instances')
91 symbols = cls._xsymbols(polyhedra)
92 islset = cls._toislset(polyhedra, symbols)
93 return cls._fromislset(islset, symbols)
94
95 @classmethod
96 def _xsymbols(cls, iterator):
97 """
98 Return the ordered tuple of symbols present in iterator.
99 """
100 symbols = set()
101 for item in iterator:
102 symbols.update(item.symbols)
103 return tuple(sorted(symbols, key=Symbol.sortkey))
104
105 @property
106 def polyhedra(self):
107 """
108 The tuple of polyhedra present in the domain.
109 """
110 return self._polyhedra
111
112 @property
113 def symbols(self):
114 """
115 The tuple of symbols present in the domain equations, sorted according
116 to Symbol.sortkey().
117 """
118 return self._symbols
119
120 @property
121 def dimension(self):
122 """
123 The dimension of the domain, i.e. the number of symbols present in it.
124 """
125 return self._dimension
126
127 def isempty(self):
128 """
129 Return True if the domain is empty, that is, equal to Empty.
130 """
131 islset = self._toislset(self.polyhedra, self.symbols)
132 empty = bool(libisl.isl_set_is_empty(islset))
133 libisl.isl_set_free(islset)
134 return empty
135
136 def __bool__(self):
137 """
138 Return True if the domain is non-empty.
139 """
140 return not self.isempty()
141
142 def isuniverse(self):
143 """
144 Return True if the domain is universal, that is, equal to Universe.
145 """
146 islset = self._toislset(self.polyhedra, self.symbols)
147 universe = bool(libisl.isl_set_plain_is_universe(islset))
148 libisl.isl_set_free(islset)
149 return universe
150
151 def isbounded(self):
152 """
153 Return True if the domain is bounded.
154 """
155 islset = self._toislset(self.polyhedra, self.symbols)
156 bounded = bool(libisl.isl_set_is_bounded(islset))
157 libisl.isl_set_free(islset)
158 return bounded
159
160 def __eq__(self, other):
161 """
162 Return True if two domains are equal.
163 """
164 symbols = self._xsymbols([self, other])
165 islset1 = self._toislset(self.polyhedra, symbols)
166 islset2 = other._toislset(other.polyhedra, symbols)
167 equal = bool(libisl.isl_set_is_equal(islset1, islset2))
168 libisl.isl_set_free(islset1)
169 libisl.isl_set_free(islset2)
170 return equal
171
172 def isdisjoint(self, other):
173 """
174 Return True if two domains have a null intersection.
175 """
176 symbols = self._xsymbols([self, other])
177 islset1 = self._toislset(self.polyhedra, symbols)
178 islset2 = self._toislset(other.polyhedra, symbols)
179 equal = bool(libisl.isl_set_is_disjoint(islset1, islset2))
180 libisl.isl_set_free(islset1)
181 libisl.isl_set_free(islset2)
182 return equal
183
184 def issubset(self, other):
185 """
186 Report whether another domain contains the domain.
187 """
188 symbols = self._xsymbols([self, other])
189 islset1 = self._toislset(self.polyhedra, symbols)
190 islset2 = self._toislset(other.polyhedra, symbols)
191 equal = bool(libisl.isl_set_is_subset(islset1, islset2))
192 libisl.isl_set_free(islset1)
193 libisl.isl_set_free(islset2)
194 return equal
195
196 def __le__(self, other):
197 return self.issubset(other)
198 __le__.__doc__ = issubset.__doc__
199
200 def __lt__(self, other):
201 """
202 Report whether another domain is contained within the domain.
203 """
204 symbols = self._xsymbols([self, other])
205 islset1 = self._toislset(self.polyhedra, symbols)
206 islset2 = self._toislset(other.polyhedra, symbols)
207 equal = bool(libisl.isl_set_is_strict_subset(islset1, islset2))
208 libisl.isl_set_free(islset1)
209 libisl.isl_set_free(islset2)
210 return equal
211
212 def complement(self):
213 """
214 Return the complementary domain of the domain.
215 """
216 islset = self._toislset(self.polyhedra, self.symbols)
217 islset = libisl.isl_set_complement(islset)
218 return self._fromislset(islset, self.symbols)
219
220 def __invert__(self):
221 return self.complement()
222 __invert__.__doc__ = complement.__doc__
223
224 def make_disjoint(self):
225 """
226 Return an equivalent domain, whose polyhedra are disjoint.
227 """
228 islset = self._toislset(self.polyhedra, self.symbols)
229 islset = libisl.isl_set_make_disjoint(mainctx, islset)
230 return self._fromislset(islset, self.symbols)
231
232 def coalesce(self):
233 """
234 Simplify the representation of the domain by trying to combine pairs of
235 polyhedra into a single polyhedron, and return the resulting domain.
236 """
237 islset = self._toislset(self.polyhedra, self.symbols)
238 islset = libisl.isl_set_coalesce(islset)
239 return self._fromislset(islset, self.symbols)
240
241 def detect_equalities(self):
242 """
243 Simplify the representation of the domain by detecting implicit
244 equalities, and return the resulting domain.
245 """
246 islset = self._toislset(self.polyhedra, self.symbols)
247 islset = libisl.isl_set_detect_equalities(islset)
248 return self._fromislset(islset, self.symbols)
249
250 def remove_redundancies(self):
251 """
252 Remove redundant constraints in the domain, and return the resulting
253 domain.
254 """
255 islset = self._toislset(self.polyhedra, self.symbols)
256 islset = libisl.isl_set_remove_redundancies(islset)
257 return self._fromislset(islset, self.symbols)
258
259 def aspolyhedron(self):
260 from .polyhedra import Polyhedron
261 islset = self._toislset(self.polyhedra, self.symbols)
262 islbset = libisl.isl_set_polyhedral_hull(islset)
263 return Polyhedron._fromislbasicset(islbset, self.symbols)
264
265 def asdomain(self):
266 return self
267
268 def project(self, symbols):
269 """
270 Project out the sequence of symbols given in arguments, and return the
271 resulting domain.
272 """
273 islset = self._toislset(self.polyhedra, self.symbols)
274 n = 0
275 for index, symbol in reversed(list(enumerate(self.symbols))):
276 if symbol in symbols:
277 n += 1
278 elif n > 0:
279 islset = libisl.isl_set_project_out(islset,
280 libisl.isl_dim_set, index + 1, n)
281 n = 0
282 if n > 0:
283 islset = libisl.isl_set_project_out(islset, libisl.isl_dim_set, 0, n)
284 symbols = [symbol for symbol in self.symbols if symbol not in symbols]
285 return Domain._fromislset(islset, symbols)
286
287 def sample(self):
288 """
289 Return a sample of the domain, as an integer instance of Point. If the
290 domain is empty, a ValueError exception is raised.
291 """
292 islset = self._toislset(self.polyhedra, self.symbols)
293 islpoint = libisl.isl_set_sample_point(islset)
294 if bool(libisl.isl_point_is_void(islpoint)):
295 libisl.isl_point_free(islpoint)
296 raise ValueError('domain must be non-empty')
297 point = {}
298 for index, symbol in enumerate(self.symbols):
299 coordinate = libisl.isl_point_get_coordinate_val(islpoint,
300 libisl.isl_dim_set, index)
301 coordinate = islhelper.isl_val_to_int(coordinate)
302 point[symbol] = coordinate
303 libisl.isl_point_free(islpoint)
304 return point
305
306 def intersection(self, *others):
307 """
308 Return the intersection of two or more domains as a new domain. As an
309 alternative, function And() can be used.
310 """
311 if len(others) == 0:
312 return self
313 symbols = self._xsymbols((self,) + others)
314 islset1 = self._toislset(self.polyhedra, symbols)
315 for other in others:
316 islset2 = other._toislset(other.polyhedra, symbols)
317 islset1 = libisl.isl_set_intersect(islset1, islset2)
318 return self._fromislset(islset1, symbols)
319
320 def __and__(self, other):
321 return self.intersection(other)
322 __and__.__doc__ = intersection.__doc__
323
324 def union(self, *others):
325 """
326 Return the union of two or more domains as a new domain. As an
327 alternative, function Or() can be used.
328 """
329 if len(others) == 0:
330 return self
331 symbols = self._xsymbols((self,) + others)
332 islset1 = self._toislset(self.polyhedra, symbols)
333 for other in others:
334 islset2 = other._toislset(other.polyhedra, symbols)
335 islset1 = libisl.isl_set_union(islset1, islset2)
336 return self._fromislset(islset1, symbols)
337
338 def __or__(self, other):
339 return self.union(other)
340 __or__.__doc__ = union.__doc__
341
342 def __add__(self, other):
343 return self.union(other)
344 __add__.__doc__ = union.__doc__
345
346 def difference(self, other):
347 """
348 Return the difference of two domains as a new domain.
349 """
350 symbols = self._xsymbols([self, other])
351 islset1 = self._toislset(self.polyhedra, symbols)
352 islset2 = other._toislset(other.polyhedra, symbols)
353 islset = libisl.isl_set_subtract(islset1, islset2)
354 return self._fromislset(islset, symbols)
355
356 def __sub__(self, other):
357 return self.difference(other)
358 __sub__.__doc__ = difference.__doc__
359
360 def lexmin(self):
361 """
362 Return the lexicographic minimum of the elements in the domain.
363 """
364 islset = self._toislset(self.polyhedra, self.symbols)
365 islset = libisl.isl_set_lexmin(islset)
366 return self._fromislset(islset, self.symbols)
367
368 def lexmax(self):
369 """
370 Return the lexicographic maximum of the elements in the domain.
371 """
372 islset = self._toislset(self.polyhedra, self.symbols)
373 islset = libisl.isl_set_lexmax(islset)
374 return self._fromislset(islset, self.symbols)
375
376 _RE_COORDINATE = re.compile(r'\((?P<num>\-?\d+)\)(/(?P<den>\d+))?')
377
378 def vertices(self):
379 """
380 Return the vertices of the domain, as a list of rational instances of
381 Point.
382 """
383 from .polyhedra import Polyhedron
384 if not self.isbounded():
385 raise ValueError('domain must be bounded')
386 islbset = self._toislbasicset(self.equalities, self.inequalities,
387 self.symbols)
388 vertices = libisl.isl_basic_set_compute_vertices(islbset);
389 vertices = islhelper.isl_vertices_vertices(vertices)
390 points = []
391 for vertex in vertices:
392 expr = libisl.isl_vertex_get_expr(vertex)
393 coordinates = []
394 if islhelper.isl_version < '0.13':
395 constraints = islhelper.isl_basic_set_constraints(expr)
396 for constraint in constraints:
397 constant = libisl.isl_constraint_get_constant_val(constraint)
398 constant = islhelper.isl_val_to_int(constant)
399 for index, symbol in enumerate(self.symbols):
400 coefficient = libisl.isl_constraint_get_coefficient_val(constraint,
401 libisl.isl_dim_set, index)
402 coefficient = islhelper.isl_val_to_int(coefficient)
403 if coefficient != 0:
404 coordinate = -Fraction(constant, coefficient)
405 coordinates.append((symbol, coordinate))
406 else:
407 string = islhelper.isl_multi_aff_to_str(expr)
408 matches = self._RE_COORDINATE.finditer(string)
409 for symbol, match in zip(self.symbols, matches):
410 numerator = int(match.group('num'))
411 denominator = match.group('den')
412 denominator = 1 if denominator is None else int(denominator)
413 coordinate = Fraction(numerator, denominator)
414 coordinates.append((symbol, coordinate))
415 points.append(Point(coordinates))
416 return points
417
418 def points(self):
419 """
420 Return the integer points of a bounded domain, as a list of integer
421 instances of Point. If the domain is not bounded, a ValueError exception
422 is raised.
423 """
424 if not self.isbounded():
425 raise ValueError('domain must be bounded')
426 from .polyhedra import Universe, Eq
427 islset = self._toislset(self.polyhedra, self.symbols)
428 islpoints = islhelper.isl_set_points(islset)
429 points = []
430 for islpoint in islpoints:
431 coordinates = {}
432 for index, symbol in enumerate(self.symbols):
433 coordinate = libisl.isl_point_get_coordinate_val(islpoint,
434 libisl.isl_dim_set, index)
435 coordinate = islhelper.isl_val_to_int(coordinate)
436 coordinates[symbol] = coordinate
437 points.append(Point(coordinates))
438 return points
439
440 def __contains__(self, point):
441 """
442 Return True if the point if contained within the domain.
443 """
444 for polyhedron in self.polyhedra:
445 if point in polyhedron:
446 return True
447 return False
448
449 @classmethod
450 def _polygon_inner_point(cls, points):
451 symbols = points[0].symbols
452 coordinates = {symbol: 0 for symbol in symbols}
453 for point in points:
454 for symbol, coordinate in point.coordinates():
455 coordinates[symbol] += coordinate
456 for symbol in symbols:
457 coordinates[symbol] /= len(points)
458 return Point(coordinates)
459
460 @classmethod
461 def _sort_polygon_2d(cls, points):
462 if len(points) <= 3:
463 return points
464 o = cls._polygon_inner_point(points)
465 angles = {}
466 for m in points:
467 om = Vector(o, m)
468 dx, dy = (coordinate for symbol, coordinate in om.coordinates())
469 angle = math.atan2(dy, dx)
470 angles[m] = angle
471 return sorted(points, key=angles.get)
472
473 @classmethod
474 def _sort_polygon_3d(cls, points):
475 if len(points) <= 3:
476 return points
477 o = cls._polygon_inner_point(points)
478 a = points[0]
479 oa = Vector(o, a)
480 norm_oa = oa.norm()
481 for b in points[1:]:
482 ob = Vector(o, b)
483 u = oa.cross(ob)
484 if not u.isnull():
485 u = u.asunit()
486 break
487 else:
488 raise ValueError('degenerate polygon')
489 angles = {a: 0.}
490 for m in points[1:]:
491 om = Vector(o, m)
492 normprod = norm_oa * om.norm()
493 cosinus = max(oa.dot(om) / normprod, -1.)
494 sinus = u.dot(oa.cross(om)) / normprod
495 angle = math.acos(cosinus)
496 angle = math.copysign(angle, sinus)
497 angles[m] = angle
498 return sorted(points, key=angles.get)
499
500 def faces(self):
501 """
502 Return the list of faces of a bounded domain. Each face is represented
503 by a list of vertices, in the form of rational instances of Point. If
504 the domain is not bounded, a ValueError exception is raised.
505 """
506 faces = []
507 for polyhedron in self.polyhedra:
508 vertices = polyhedron.vertices()
509 for constraint in polyhedron.constraints:
510 face = []
511 for vertex in vertices:
512 if constraint.subs(vertex.coordinates()) == 0:
513 face.append(vertex)
514 if len(face) >= 3:
515 faces.append(face)
516 return faces
517
518 def _plot_2d(self, plot=None, **kwargs):
519 import matplotlib.pyplot as plt
520 from matplotlib.patches import Polygon
521 if plot is None:
522 fig = plt.figure()
523 plot = fig.add_subplot(1, 1, 1)
524 xmin, xmax = plot.get_xlim()
525 ymin, ymax = plot.get_ylim()
526 for polyhedron in self.polyhedra:
527 vertices = polyhedron._sort_polygon_2d(polyhedron.vertices())
528 xys = [tuple(vertex.values()) for vertex in vertices]
529 xs, ys = zip(*xys)
530 xmin, xmax = min(xmin, float(min(xs))), max(xmax, float(max(xs)))
531 ymin, ymax = min(ymin, float(min(ys))), max(ymax, float(max(ys)))
532 plot.add_patch(Polygon(xys, closed=True, **kwargs))
533 plot.set_xlim(xmin, xmax)
534 plot.set_ylim(ymin, ymax)
535 return plot
536
537 def _plot_3d(self, plot=None, **kwargs):
538 import matplotlib.pyplot as plt
539 from mpl_toolkits.mplot3d import Axes3D
540 from mpl_toolkits.mplot3d.art3d import Poly3DCollection
541 if plot is None:
542 fig = plt.figure()
543 axes = Axes3D(fig)
544 else:
545 axes = plot
546 xmin, xmax = axes.get_xlim()
547 ymin, ymax = axes.get_ylim()
548 zmin, zmax = axes.get_zlim()
549 poly_xyzs = []
550 for vertices in self.faces():
551 vertices = self._sort_polygon_3d(vertices)
552 vertices.append(vertices[0])
553 face_xyzs = [tuple(vertex.values()) for vertex in vertices]
554 xs, ys, zs = zip(*face_xyzs)
555 xmin, xmax = min(xmin, float(min(xs))), max(xmax, float(max(xs)))
556 ymin, ymax = min(ymin, float(min(ys))), max(ymax, float(max(ys)))
557 zmin, zmax = min(zmin, float(min(zs))), max(zmax, float(max(zs)))
558 poly_xyzs.append(face_xyzs)
559 collection = Poly3DCollection(poly_xyzs, **kwargs)
560 axes.add_collection3d(collection)
561 axes.set_xlim(xmin, xmax)
562 axes.set_ylim(ymin, ymax)
563 axes.set_zlim(zmin, zmax)
564 return axes
565
566 def plot(self, plot=None, **kwargs):
567 """
568 Plot a 2D or 3D domain using matplotlib. Draw it to the current plot
569 object if present, otherwise create a new one. options are keyword
570 arguments passed to the matplotlib drawing functions, they can be used
571 to set the drawing color for example. Raise ValueError is the domain is
572 not 2D or 3D.
573 """
574 if not self.isbounded():
575 raise ValueError('domain must be bounded')
576 elif self.dimension == 2:
577 return self._plot_2d(plot=plot, **kwargs)
578 elif self.dimension == 3:
579 return self._plot_3d(plot=plot, **kwargs)
580 else:
581 raise ValueError('polyhedron must be 2 or 3-dimensional')
582
583 def subs(self, symbol, expression=None):
584 """
585 Substitute the given symbol by an expression in the domain constraints.
586 To perform multiple substitutions at once, pass a sequence or a
587 dictionary of (old, new) pairs to subs. The syntax of this function is
588 similar to LinExpr.subs().
589 """
590 polyhedra = [polyhedron.subs(symbol, expression)
591 for polyhedron in self.polyhedra]
592 return Domain(*polyhedra)
593
594 @classmethod
595 def _fromislset(cls, islset, symbols):
596 from .polyhedra import Polyhedron
597 islset = libisl.isl_set_remove_divs(islset)
598 islbsets = islhelper.isl_set_basic_sets(islset)
599 libisl.isl_set_free(islset)
600 polyhedra = []
601 for islbset in islbsets:
602 polyhedron = Polyhedron._fromislbasicset(islbset, symbols)
603 polyhedra.append(polyhedron)
604 if len(polyhedra) == 0:
605 from .polyhedra import Empty
606 return Empty
607 elif len(polyhedra) == 1:
608 return polyhedra[0]
609 else:
610 self = object().__new__(Domain)
611 self._polyhedra = tuple(polyhedra)
612 self._symbols = cls._xsymbols(polyhedra)
613 self._dimension = len(self._symbols)
614 return self
615
616 @classmethod
617 def _toislset(cls, polyhedra, symbols):
618 polyhedron = polyhedra[0]
619 islbset = polyhedron._toislbasicset(polyhedron.equalities,
620 polyhedron.inequalities, symbols)
621 islset1 = libisl.isl_set_from_basic_set(islbset)
622 for polyhedron in polyhedra[1:]:
623 islbset = polyhedron._toislbasicset(polyhedron.equalities,
624 polyhedron.inequalities, symbols)
625 islset2 = libisl.isl_set_from_basic_set(islbset)
626 islset1 = libisl.isl_set_union(islset1, islset2)
627 return islset1
628
629 @classmethod
630 def _fromast(cls, node):
631 from .polyhedra import Polyhedron
632 if isinstance(node, ast.Module) and len(node.body) == 1:
633 return cls._fromast(node.body[0])
634 elif isinstance(node, ast.Expr):
635 return cls._fromast(node.value)
636 elif isinstance(node, ast.UnaryOp):
637 domain = cls._fromast(node.operand)
638 if isinstance(node.operand, ast.invert):
639 return Not(domain)
640 elif isinstance(node, ast.BinOp):
641 domain1 = cls._fromast(node.left)
642 domain2 = cls._fromast(node.right)
643 if isinstance(node.op, ast.BitAnd):
644 return And(domain1, domain2)
645 elif isinstance(node.op, ast.BitOr):
646 return Or(domain1, domain2)
647 elif isinstance(node, ast.Compare):
648 equalities = []
649 inequalities = []
650 left = LinExpr._fromast(node.left)
651 for i in range(len(node.ops)):
652 op = node.ops[i]
653 right = LinExpr._fromast(node.comparators[i])
654 if isinstance(op, ast.Lt):
655 inequalities.append(right - left - 1)
656 elif isinstance(op, ast.LtE):
657 inequalities.append(right - left)
658 elif isinstance(op, ast.Eq):
659 equalities.append(left - right)
660 elif isinstance(op, ast.GtE):
661 inequalities.append(left - right)
662 elif isinstance(op, ast.Gt):
663 inequalities.append(left - right - 1)
664 else:
665 break
666 left = right
667 else:
668 return Polyhedron(equalities, inequalities)
669 raise SyntaxError('invalid syntax')
670
671 _RE_BRACES = re.compile(r'^\{\s*|\s*\}$')
672 _RE_EQ = re.compile(r'([^<=>])=([^<=>])')
673 _RE_AND = re.compile(r'\band\b|,|&&|/\\|∧|∩')
674 _RE_OR = re.compile(r'\bor\b|;|\|\||\\/|∨|∪')
675 _RE_NOT = re.compile(r'\bnot\b|!|¬')
676 _RE_NUM_VAR = LinExpr._RE_NUM_VAR
677 _RE_OPERATORS = re.compile(r'(&|\||~)')
678
679 @classmethod
680 def fromstring(cls, string):
681 """
682 Create a domain from a string. Raise SyntaxError if the string is not
683 properly formatted.
684 """
685 # remove curly brackets
686 string = cls._RE_BRACES.sub(r'', string)
687 # replace '=' by '=='
688 string = cls._RE_EQ.sub(r'\1==\2', string)
689 # replace 'and', 'or', 'not'
690 string = cls._RE_AND.sub(r' & ', string)
691 string = cls._RE_OR.sub(r' | ', string)
692 string = cls._RE_NOT.sub(r' ~', string)
693 # add implicit multiplication operators, e.g. '5x' -> '5*x'
694 string = cls._RE_NUM_VAR.sub(r'\1*\2', string)
695 # add parentheses to force precedence
696 tokens = cls._RE_OPERATORS.split(string)
697 for i, token in enumerate(tokens):
698 if i % 2 == 0:
699 token = '({})'.format(token)
700 tokens[i] = token
701 string = ''.join(tokens)
702 tree = ast.parse(string, 'eval')
703 return cls._fromast(tree)
704
705 def __repr__(self):
706 assert len(self.polyhedra) >= 2
707 strings = [repr(polyhedron) for polyhedron in self.polyhedra]
708 return 'Or({})'.format(', '.join(strings))
709
710 def _repr_latex_(self):
711 strings = []
712 for polyhedron in self.polyhedra:
713 strings.append('({})'.format(polyhedron._repr_latex_().strip('$')))
714 return '${}$'.format(' \\vee '.join(strings))
715
716 @classmethod
717 def fromsympy(cls, expr):
718 """
719 Create a domain from a sympy expression.
720 """
721 import sympy
722 from .polyhedra import Lt, Le, Eq, Ne, Ge, Gt
723 funcmap = {
724 sympy.And: And, sympy.Or: Or, sympy.Not: Not,
725 sympy.Lt: Lt, sympy.Le: Le,
726 sympy.Eq: Eq, sympy.Ne: Ne,
727 sympy.Ge: Ge, sympy.Gt: Gt,
728 }
729 if expr.func in funcmap:
730 args = [Domain.fromsympy(arg) for arg in expr.args]
731 return funcmap[expr.func](*args)
732 elif isinstance(expr, sympy.Expr):
733 return LinExpr.fromsympy(expr)
734 raise ValueError('non-domain expression: {!r}'.format(expr))
735
736 def tosympy(self):
737 """
738 Convert the domain to a sympy expression.
739 """
740 import sympy
741 polyhedra = [polyhedron.tosympy() for polyhedron in polyhedra]
742 return sympy.Or(*polyhedra)
743
744
745 def And(*domains):
746 """
747 Create the intersection domain of the domains given in arguments.
748 """
749 if len(domains) == 0:
750 from .polyhedra import Universe
751 return Universe
752 else:
753 return domains[0].intersection(*domains[1:])
754 And.__doc__ = Domain.intersection.__doc__
755
756 def Or(*domains):
757 """
758 Create the union domain of the domains given in arguments.
759 """
760 if len(domains) == 0:
761 from .polyhedra import Empty
762 return Empty
763 else:
764 return domains[0].union(*domains[1:])
765 Or.__doc__ = Domain.union.__doc__
766
767 def Not(domain):
768 """
769 Create the complementary domain of the domain given in argument.
770 """
771 return ~domain
772 Not.__doc__ = Domain.complement.__doc__