da0ea97d945510b44e0cfa6310a406de8e9123eb
[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 if isinstance(other, Domain):
165 symbols = self._xsymbols([self, other])
166 islset1 = self._toislset(self.polyhedra, symbols)
167 islset2 = other._toislset(other.polyhedra, symbols)
168 equal = bool(libisl.isl_set_is_equal(islset1, islset2))
169 libisl.isl_set_free(islset1)
170 libisl.isl_set_free(islset2)
171 return equal
172 return NotImplemented
173
174 def isdisjoint(self, other):
175 """
176 Return True if two domains have a null intersection.
177 """
178 if not isinstance(other, Domain):
179 raise TypeError('other must be a Domain instance')
180 symbols = self._xsymbols([self, other])
181 islset1 = self._toislset(self.polyhedra, symbols)
182 islset2 = self._toislset(other.polyhedra, symbols)
183 equal = bool(libisl.isl_set_is_disjoint(islset1, islset2))
184 libisl.isl_set_free(islset1)
185 libisl.isl_set_free(islset2)
186 return equal
187
188 def issubset(self, other):
189 """
190 Report whether another domain contains the domain.
191 """
192 return self < other
193
194 def __le__(self, other):
195 if isinstance(other, Domain):
196 symbols = self._xsymbols([self, other])
197 islset1 = self._toislset(self.polyhedra, symbols)
198 islset2 = self._toislset(other.polyhedra, symbols)
199 equal = bool(libisl.isl_set_is_subset(islset1, islset2))
200 libisl.isl_set_free(islset1)
201 libisl.isl_set_free(islset2)
202 return equal
203 return NotImplemented
204 __le__.__doc__ = issubset.__doc__
205
206 def __lt__(self, other):
207 """
208 Report whether another domain is contained within the domain.
209 """
210 if isinstance(other, Domain):
211 symbols = self._xsymbols([self, other])
212 islset1 = self._toislset(self.polyhedra, symbols)
213 islset2 = self._toislset(other.polyhedra, symbols)
214 equal = bool(libisl.isl_set_is_strict_subset(islset1, islset2))
215 libisl.isl_set_free(islset1)
216 libisl.isl_set_free(islset2)
217 return equal
218 return NotImplemented
219
220 def complement(self):
221 """
222 Return the complementary domain of the domain.
223 """
224 islset = self._toislset(self.polyhedra, self.symbols)
225 islset = libisl.isl_set_complement(islset)
226 return self._fromislset(islset, self.symbols)
227
228 def __invert__(self):
229 return self.complement()
230 __invert__.__doc__ = complement.__doc__
231
232 def make_disjoint(self):
233 """
234 Return an equivalent domain, whose polyhedra are disjoint.
235 """
236 islset = self._toislset(self.polyhedra, self.symbols)
237 islset = libisl.isl_set_make_disjoint(mainctx, islset)
238 return self._fromislset(islset, self.symbols)
239
240 def coalesce(self):
241 """
242 Simplify the representation of the domain by trying to combine pairs of
243 polyhedra into a single polyhedron, and return the resulting domain.
244 """
245 islset = self._toislset(self.polyhedra, self.symbols)
246 islset = libisl.isl_set_coalesce(islset)
247 return self._fromislset(islset, self.symbols)
248
249 def detect_equalities(self):
250 """
251 Simplify the representation of the domain by detecting implicit
252 equalities, and return the resulting domain.
253 """
254 islset = self._toislset(self.polyhedra, self.symbols)
255 islset = libisl.isl_set_detect_equalities(islset)
256 return self._fromislset(islset, self.symbols)
257
258 def remove_redundancies(self):
259 """
260 Remove redundant constraints in the domain, and return the resulting
261 domain.
262 """
263 islset = self._toislset(self.polyhedra, self.symbols)
264 islset = libisl.isl_set_remove_redundancies(islset)
265 return self._fromislset(islset, self.symbols)
266
267 def aspolyhedron(self):
268 from .polyhedra import Polyhedron
269 islset = self._toislset(self.polyhedra, self.symbols)
270 islbset = libisl.isl_set_polyhedral_hull(islset)
271 return Polyhedron._fromislbasicset(islbset, self.symbols)
272
273 def asdomain(self):
274 return self
275
276 def project(self, symbols):
277 """
278 Project out the sequence of symbols given in arguments, and return the
279 resulting domain.
280 """
281 symbols = list(symbols)
282 for symbol in symbols:
283 if not isinstance(symbol, Symbol):
284 raise TypeError('symbols must be Symbol instances')
285 islset = self._toislset(self.polyhedra, self.symbols)
286 n = 0
287 for index, symbol in reversed(list(enumerate(self.symbols))):
288 if symbol in symbols:
289 n += 1
290 elif n > 0:
291 islset = libisl.isl_set_project_out(islset,
292 libisl.isl_dim_set, index + 1, n)
293 n = 0
294 if n > 0:
295 islset = libisl.isl_set_project_out(islset, libisl.isl_dim_set, 0, n)
296 symbols = [symbol for symbol in self.symbols if symbol not in symbols]
297 return Domain._fromislset(islset, symbols)
298
299 def sample(self):
300 """
301 Return a sample of the domain, as an integer instance of Point. If the
302 domain is empty, a ValueError exception is raised.
303 """
304 islset = self._toislset(self.polyhedra, self.symbols)
305 islpoint = libisl.isl_set_sample_point(islset)
306 if bool(libisl.isl_point_is_void(islpoint)):
307 libisl.isl_point_free(islpoint)
308 raise ValueError('domain must be non-empty')
309 point = {}
310 for index, symbol in enumerate(self.symbols):
311 coordinate = libisl.isl_point_get_coordinate_val(islpoint,
312 libisl.isl_dim_set, index)
313 coordinate = islhelper.isl_val_to_int(coordinate)
314 point[symbol] = coordinate
315 libisl.isl_point_free(islpoint)
316 return point
317
318 def intersection(self, *others):
319 """
320 Return the intersection of two or more domains as a new domain. As an
321 alternative, function And() can be used.
322 """
323 result = self
324 for other in others:
325 result &= other
326 return result
327
328 def __and__(self, other):
329 if isinstance(other, Domain):
330 symbols = self._xsymbols([self, other])
331 islset1 = self._toislset(self.polyhedra, symbols)
332 islset2 = other._toislset(other.polyhedra, symbols)
333 islset = libisl.isl_set_intersect(islset1, islset2)
334 return self._fromislset(islset, symbols)
335 return NotImplemented
336 __and__.__doc__ = intersection.__doc__
337
338 def union(self, *others):
339 """
340 Return the union of two or more domains as a new domain. As an
341 alternative, function Or() can be used.
342 """
343 if len(others) == 0:
344 return self
345 symbols = self._xsymbols((self,) + others)
346 islset1 = self._toislset(self.polyhedra, symbols)
347 for other in others:
348 islset2 = other._toislset(other.polyhedra, symbols)
349 islset1 = libisl.isl_set_union(islset1, islset2)
350 return self._fromislset(islset1, symbols)
351
352 def __or__(self, other):
353 return self.union(other)
354 __or__.__doc__ = union.__doc__
355
356 def __add__(self, other):
357 return self.union(other)
358 __add__.__doc__ = union.__doc__
359
360 def difference(self, other):
361 """
362 Return the difference of two domains as a new domain.
363 """
364 symbols = self._xsymbols([self, other])
365 islset1 = self._toislset(self.polyhedra, symbols)
366 islset2 = other._toislset(other.polyhedra, symbols)
367 islset = libisl.isl_set_subtract(islset1, islset2)
368 return self._fromislset(islset, symbols)
369
370 def __sub__(self, other):
371 return self.difference(other)
372 __sub__.__doc__ = difference.__doc__
373
374 def lexmin(self):
375 """
376 Return the lexicographic minimum of the elements in the domain.
377 """
378 islset = self._toislset(self.polyhedra, self.symbols)
379 islset = libisl.isl_set_lexmin(islset)
380 return self._fromislset(islset, self.symbols)
381
382 def lexmax(self):
383 """
384 Return the lexicographic maximum of the elements in the domain.
385 """
386 islset = self._toislset(self.polyhedra, self.symbols)
387 islset = libisl.isl_set_lexmax(islset)
388 return self._fromislset(islset, self.symbols)
389
390 _RE_COORDINATE = re.compile(r'\((?P<num>\-?\d+)\)(/(?P<den>\d+))?')
391
392 def vertices(self):
393 """
394 Return the vertices of the domain, as a list of rational instances of
395 Point.
396 """
397 from .polyhedra import Polyhedron
398 if not self.isbounded():
399 raise ValueError('domain must be bounded')
400 islbset = self._toislbasicset(self.equalities, self.inequalities,
401 self.symbols)
402 vertices = libisl.isl_basic_set_compute_vertices(islbset);
403 vertices = islhelper.isl_vertices_vertices(vertices)
404 points = []
405 for vertex in vertices:
406 expr = libisl.isl_vertex_get_expr(vertex)
407 coordinates = []
408 if islhelper.isl_version < '0.13':
409 constraints = islhelper.isl_basic_set_constraints(expr)
410 for constraint in constraints:
411 constant = libisl.isl_constraint_get_constant_val(constraint)
412 constant = islhelper.isl_val_to_int(constant)
413 for index, symbol in enumerate(self.symbols):
414 coefficient = libisl.isl_constraint_get_coefficient_val(constraint,
415 libisl.isl_dim_set, index)
416 coefficient = islhelper.isl_val_to_int(coefficient)
417 if coefficient != 0:
418 coordinate = -Fraction(constant, coefficient)
419 coordinates.append((symbol, coordinate))
420 else:
421 string = islhelper.isl_multi_aff_to_str(expr)
422 matches = self._RE_COORDINATE.finditer(string)
423 for symbol, match in zip(self.symbols, matches):
424 numerator = int(match.group('num'))
425 denominator = match.group('den')
426 denominator = 1 if denominator is None else int(denominator)
427 coordinate = Fraction(numerator, denominator)
428 coordinates.append((symbol, coordinate))
429 points.append(Point(coordinates))
430 return points
431
432 def points(self):
433 """
434 Return the integer points of a bounded domain, as a list of integer
435 instances of Point. If the domain is not bounded, a ValueError exception
436 is raised.
437 """
438 if not self.isbounded():
439 raise ValueError('domain must be bounded')
440 from .polyhedra import Universe, Eq
441 islset = self._toislset(self.polyhedra, self.symbols)
442 islpoints = islhelper.isl_set_points(islset)
443 points = []
444 for islpoint in islpoints:
445 coordinates = {}
446 for index, symbol in enumerate(self.symbols):
447 coordinate = libisl.isl_point_get_coordinate_val(islpoint,
448 libisl.isl_dim_set, index)
449 coordinate = islhelper.isl_val_to_int(coordinate)
450 coordinates[symbol] = coordinate
451 points.append(Point(coordinates))
452 return points
453
454 def __contains__(self, point):
455 """
456 Return True if the point if contained within the domain.
457 """
458 for polyhedron in self.polyhedra:
459 if point in polyhedron:
460 return True
461 return False
462
463 @classmethod
464 def _polygon_inner_point(cls, points):
465 symbols = points[0].symbols
466 coordinates = {symbol: 0 for symbol in symbols}
467 for point in points:
468 for symbol, coordinate in point.coordinates():
469 coordinates[symbol] += coordinate
470 for symbol in symbols:
471 coordinates[symbol] /= len(points)
472 return Point(coordinates)
473
474 @classmethod
475 def _sort_polygon_2d(cls, points):
476 if len(points) <= 3:
477 return points
478 o = cls._polygon_inner_point(points)
479 angles = {}
480 for m in points:
481 om = Vector(o, m)
482 dx, dy = (coordinate for symbol, coordinate in om.coordinates())
483 angle = math.atan2(dy, dx)
484 angles[m] = angle
485 return sorted(points, key=angles.get)
486
487 @classmethod
488 def _sort_polygon_3d(cls, points):
489 if len(points) <= 3:
490 return points
491 o = cls._polygon_inner_point(points)
492 a = points[0]
493 oa = Vector(o, a)
494 norm_oa = oa.norm()
495 for b in points[1:]:
496 ob = Vector(o, b)
497 u = oa.cross(ob)
498 if not u.isnull():
499 u = u.asunit()
500 break
501 else:
502 raise ValueError('degenerate polygon')
503 angles = {a: 0.}
504 for m in points[1:]:
505 om = Vector(o, m)
506 normprod = norm_oa * om.norm()
507 cosinus = max(oa.dot(om) / normprod, -1.)
508 sinus = u.dot(oa.cross(om)) / normprod
509 angle = math.acos(cosinus)
510 angle = math.copysign(angle, sinus)
511 angles[m] = angle
512 return sorted(points, key=angles.get)
513
514 def faces(self):
515 """
516 Return the list of faces of a bounded domain. Each face is represented
517 by a list of vertices, in the form of rational instances of Point. If
518 the domain is not bounded, a ValueError exception is raised.
519 """
520 faces = []
521 for polyhedron in self.polyhedra:
522 vertices = polyhedron.vertices()
523 for constraint in polyhedron.constraints:
524 face = []
525 for vertex in vertices:
526 if constraint.subs(vertex.coordinates()) == 0:
527 face.append(vertex)
528 if len(face) >= 3:
529 faces.append(face)
530 return faces
531
532 def _plot_2d(self, plot=None, **kwargs):
533 import matplotlib.pyplot as plt
534 from matplotlib.patches import Polygon
535 if plot is None:
536 fig = plt.figure()
537 plot = fig.add_subplot(1, 1, 1)
538 xmin, xmax = plot.get_xlim()
539 ymin, ymax = plot.get_ylim()
540 for polyhedron in self.polyhedra:
541 vertices = polyhedron._sort_polygon_2d(polyhedron.vertices())
542 xys = [tuple(vertex.values()) for vertex in vertices]
543 xs, ys = zip(*xys)
544 xmin, xmax = min(xmin, float(min(xs))), max(xmax, float(max(xs)))
545 ymin, ymax = min(ymin, float(min(ys))), max(ymax, float(max(ys)))
546 plot.add_patch(Polygon(xys, closed=True, **kwargs))
547 plot.set_xlim(xmin, xmax)
548 plot.set_ylim(ymin, ymax)
549 return plot
550
551 def _plot_3d(self, plot=None, **kwargs):
552 import matplotlib.pyplot as plt
553 from mpl_toolkits.mplot3d import Axes3D
554 from mpl_toolkits.mplot3d.art3d import Poly3DCollection
555 if plot is None:
556 fig = plt.figure()
557 axes = Axes3D(fig)
558 else:
559 axes = plot
560 xmin, xmax = axes.get_xlim()
561 ymin, ymax = axes.get_ylim()
562 zmin, zmax = axes.get_zlim()
563 poly_xyzs = []
564 for vertices in self.faces():
565 vertices = self._sort_polygon_3d(vertices)
566 vertices.append(vertices[0])
567 face_xyzs = [tuple(vertex.values()) for vertex in vertices]
568 xs, ys, zs = zip(*face_xyzs)
569 xmin, xmax = min(xmin, float(min(xs))), max(xmax, float(max(xs)))
570 ymin, ymax = min(ymin, float(min(ys))), max(ymax, float(max(ys)))
571 zmin, zmax = min(zmin, float(min(zs))), max(zmax, float(max(zs)))
572 poly_xyzs.append(face_xyzs)
573 collection = Poly3DCollection(poly_xyzs, **kwargs)
574 axes.add_collection3d(collection)
575 axes.set_xlim(xmin, xmax)
576 axes.set_ylim(ymin, ymax)
577 axes.set_zlim(zmin, zmax)
578 return axes
579
580 def plot(self, plot=None, **kwargs):
581 """
582 Plot a 2D or 3D domain using matplotlib. Draw it to the current plot
583 object if present, otherwise create a new one. options are keyword
584 arguments passed to the matplotlib drawing functions, they can be used
585 to set the drawing color for example. Raise ValueError is the domain is
586 not 2D or 3D.
587 """
588 if not self.isbounded():
589 raise ValueError('domain must be bounded')
590 elif self.dimension == 2:
591 return self._plot_2d(plot=plot, **kwargs)
592 elif self.dimension == 3:
593 return self._plot_3d(plot=plot, **kwargs)
594 else:
595 raise ValueError('polyhedron must be 2 or 3-dimensional')
596
597 def subs(self, symbol, expression=None):
598 """
599 Substitute the given symbol by an expression in the domain constraints.
600 To perform multiple substitutions at once, pass a sequence or a
601 dictionary of (old, new) pairs to subs. The syntax of this function is
602 similar to LinExpr.subs().
603 """
604 polyhedra = [polyhedron.subs(symbol, expression)
605 for polyhedron in self.polyhedra]
606 return Domain(*polyhedra)
607
608 @classmethod
609 def _fromislset(cls, islset, symbols):
610 from .polyhedra import Polyhedron
611 islset = libisl.isl_set_remove_divs(islset)
612 islbsets = islhelper.isl_set_basic_sets(islset)
613 libisl.isl_set_free(islset)
614 polyhedra = []
615 for islbset in islbsets:
616 polyhedron = Polyhedron._fromislbasicset(islbset, symbols)
617 polyhedra.append(polyhedron)
618 if len(polyhedra) == 0:
619 from .polyhedra import Empty
620 return Empty
621 elif len(polyhedra) == 1:
622 return polyhedra[0]
623 else:
624 self = object().__new__(Domain)
625 self._polyhedra = tuple(polyhedra)
626 self._symbols = cls._xsymbols(polyhedra)
627 self._dimension = len(self._symbols)
628 return self
629
630 @classmethod
631 def _toislset(cls, polyhedra, symbols):
632 polyhedron = polyhedra[0]
633 islbset = polyhedron._toislbasicset(polyhedron.equalities,
634 polyhedron.inequalities, symbols)
635 islset1 = libisl.isl_set_from_basic_set(islbset)
636 for polyhedron in polyhedra[1:]:
637 islbset = polyhedron._toislbasicset(polyhedron.equalities,
638 polyhedron.inequalities, symbols)
639 islset2 = libisl.isl_set_from_basic_set(islbset)
640 islset1 = libisl.isl_set_union(islset1, islset2)
641 return islset1
642
643 @classmethod
644 def _fromast(cls, node):
645 from .polyhedra import Polyhedron
646 if isinstance(node, ast.Module) and len(node.body) == 1:
647 return cls._fromast(node.body[0])
648 elif isinstance(node, ast.Expr):
649 return cls._fromast(node.value)
650 elif isinstance(node, ast.UnaryOp):
651 domain = cls._fromast(node.operand)
652 if isinstance(node.operand, ast.invert):
653 return Not(domain)
654 elif isinstance(node, ast.BinOp):
655 domain1 = cls._fromast(node.left)
656 domain2 = cls._fromast(node.right)
657 if isinstance(node.op, ast.BitAnd):
658 return And(domain1, domain2)
659 elif isinstance(node.op, ast.BitOr):
660 return Or(domain1, domain2)
661 elif isinstance(node, ast.Compare):
662 equalities = []
663 inequalities = []
664 left = LinExpr._fromast(node.left)
665 for i in range(len(node.ops)):
666 op = node.ops[i]
667 right = LinExpr._fromast(node.comparators[i])
668 if isinstance(op, ast.Lt):
669 inequalities.append(right - left - 1)
670 elif isinstance(op, ast.LtE):
671 inequalities.append(right - left)
672 elif isinstance(op, ast.Eq):
673 equalities.append(left - right)
674 elif isinstance(op, ast.GtE):
675 inequalities.append(left - right)
676 elif isinstance(op, ast.Gt):
677 inequalities.append(left - right - 1)
678 else:
679 break
680 left = right
681 else:
682 return Polyhedron(equalities, inequalities)
683 raise SyntaxError('invalid syntax')
684
685 _RE_BRACES = re.compile(r'^\{\s*|\s*\}$')
686 _RE_EQ = re.compile(r'([^<=>])=([^<=>])')
687 _RE_AND = re.compile(r'\band\b|,|&&|/\\|∧|∩')
688 _RE_OR = re.compile(r'\bor\b|;|\|\||\\/|∨|∪')
689 _RE_NOT = re.compile(r'\bnot\b|!|¬')
690 _RE_NUM_VAR = LinExpr._RE_NUM_VAR
691 _RE_OPERATORS = re.compile(r'(&|\||~)')
692
693 @classmethod
694 def fromstring(cls, string):
695 """
696 Create a domain from a string. Raise SyntaxError if the string is not
697 properly formatted.
698 """
699 # remove curly brackets
700 string = cls._RE_BRACES.sub(r'', string)
701 # replace '=' by '=='
702 string = cls._RE_EQ.sub(r'\1==\2', string)
703 # replace 'and', 'or', 'not'
704 string = cls._RE_AND.sub(r' & ', string)
705 string = cls._RE_OR.sub(r' | ', string)
706 string = cls._RE_NOT.sub(r' ~', string)
707 # add implicit multiplication operators, e.g. '5x' -> '5*x'
708 string = cls._RE_NUM_VAR.sub(r'\1*\2', string)
709 # add parentheses to force precedence
710 tokens = cls._RE_OPERATORS.split(string)
711 for i, token in enumerate(tokens):
712 if i % 2 == 0:
713 token = '({})'.format(token)
714 tokens[i] = token
715 string = ''.join(tokens)
716 tree = ast.parse(string, 'eval')
717 return cls._fromast(tree)
718
719 def __repr__(self):
720 assert len(self.polyhedra) >= 2
721 strings = [repr(polyhedron) for polyhedron in self.polyhedra]
722 return 'Or({})'.format(', '.join(strings))
723
724 def _repr_latex_(self):
725 strings = []
726 for polyhedron in self.polyhedra:
727 strings.append('({})'.format(polyhedron._repr_latex_().strip('$')))
728 return '${}$'.format(' \\vee '.join(strings))
729
730 @classmethod
731 def fromsympy(cls, expr):
732 """
733 Create a domain from a sympy expression.
734 """
735 import sympy
736 from .polyhedra import Lt, Le, Eq, Ne, Ge, Gt
737 funcmap = {
738 sympy.And: And, sympy.Or: Or, sympy.Not: Not,
739 sympy.Lt: Lt, sympy.Le: Le,
740 sympy.Eq: Eq, sympy.Ne: Ne,
741 sympy.Ge: Ge, sympy.Gt: Gt,
742 }
743 if expr.func in funcmap:
744 args = [Domain.fromsympy(arg) for arg in expr.args]
745 return funcmap[expr.func](*args)
746 elif isinstance(expr, sympy.Expr):
747 return LinExpr.fromsympy(expr)
748 raise ValueError('non-domain expression: {!r}'.format(expr))
749
750 def tosympy(self):
751 """
752 Convert the domain to a sympy expression.
753 """
754 import sympy
755 polyhedra = [polyhedron.tosympy() for polyhedron in polyhedra]
756 return sympy.Or(*polyhedra)
757
758
759 def And(*domains):
760 """
761 Create the intersection domain of the domains given in arguments.
762 """
763 if len(domains) == 0:
764 from .polyhedra import Universe
765 return Universe
766 else:
767 return domains[0].intersection(*domains[1:])
768 And.__doc__ = Domain.intersection.__doc__
769
770 def Or(*domains):
771 """
772 Create the union domain of the domains given in arguments.
773 """
774 if len(domains) == 0:
775 from .polyhedra import Empty
776 return Empty
777 else:
778 return domains[0].union(*domains[1:])
779 Or.__doc__ = Domain.union.__doc__
780
781 def Not(domain):
782 """
783 Create the complementary domain of the domain given in argument.
784 """
785 return ~domain
786 Not.__doc__ = Domain.complement.__doc__