X-Git-Url: https://svn.cri.ensmp.fr/git/linpy.git/blobdiff_plain/7b93cea1daf2889e9ee10ca9c22a1b5124404937..430dff79f35801b018bfa1af2cc293f1a230a737:/linpy/domains.py?ds=inline diff --git a/linpy/domains.py b/linpy/domains.py index 21e78db..f26b2fc 100644 --- a/linpy/domains.py +++ b/linpy/domains.py @@ -24,7 +24,7 @@ from fractions import Fraction from . import islhelper from .islhelper import mainctx, libisl -from .linexprs import Expression, Symbol, Rational +from .linexprs import LinExpr, Symbol from .geometry import GeometricObject, Point, Vector @@ -36,6 +36,13 @@ __all__ = [ @functools.total_ordering class Domain(GeometricObject): + """ + A domain is a union of polyhedra. Unlike polyhedra, domains allow exact + computation of union, subtraction and complementary operations. + + A domain with a unique polyhedron is automatically subclassed as a + Polyhedron instance. + """ __slots__ = ( '_polyhedra', @@ -44,6 +51,30 @@ class Domain(GeometricObject): ) def __new__(cls, *polyhedra): + """ + Return a domain from a sequence of polyhedra. + + >>> square1 = Polyhedron('0 <= x <= 2, 0 <= y <= 2') + >>> square2 = Polyhedron('1 <= x <= 3, 1 <= y <= 3') + >>> dom = Domain(square1, square2) + >>> dom + Or(And(x <= 2, 0 <= x, y <= 2, 0 <= y), + And(x <= 3, 1 <= x, y <= 3, 1 <= y)) + + It is also possible to build domains from polyhedra using arithmetic + operators Domain.__or__(), Domain.__invert__() or functions Or() and + Not(), using one of the following instructions: + + >>> dom = square1 | square2 + >>> dom = Or(square1, square2) + + Alternatively, a domain can be built from a string: + + >>> dom = Domain('0 <= x <= 2, 0 <= y <= 2; 1 <= x <= 3, 1 <= y <= 3') + + Finally, a domain can be built from a GeometricObject instance, calling + the GeometricObject.asdomain() method. + """ from .polyhedra import Polyhedron if len(polyhedra) == 1: argument = polyhedra[0] @@ -74,27 +105,29 @@ class Domain(GeometricObject): @property def polyhedra(self): + """ + The tuple of polyhedra present in the domain. + """ return self._polyhedra @property def symbols(self): + """ + The tuple of symbols present in the domain equations, sorted according + to Symbol.sortkey(). + """ return self._symbols @property def dimension(self): - return self._dimension - - def disjoint(self): """ - Returns this set as disjoint. + The dimension of the domain, i.e. the number of symbols present in it. """ - islset = self._toislset(self.polyhedra, self.symbols) - islset = libisl.isl_set_make_disjoint(mainctx, islset) - return self._fromislset(islset, self.symbols) + return self._dimension def isempty(self): """ - Returns true if this set is an Empty set. + Return True if the domain is empty, that is, equal to Empty. """ islset = self._toislset(self.polyhedra, self.symbols) empty = bool(libisl.isl_set_is_empty(islset)) @@ -102,11 +135,14 @@ class Domain(GeometricObject): return empty def __bool__(self): + """ + Return True if the domain is non-empty. + """ return not self.isempty() def isuniverse(self): """ - Returns true if this set is the Universe set. + Return True if the domain is universal, that is, equal to Universe. """ islset = self._toislset(self.polyhedra, self.symbols) universe = bool(libisl.isl_set_plain_is_universe(islset)) @@ -115,7 +151,7 @@ class Domain(GeometricObject): def isbounded(self): """ - Returns true if this set is bounded. + Return True if the domain is bounded. """ islset = self._toislset(self.polyhedra, self.symbols) bounded = bool(libisl.isl_set_is_bounded(islset)) @@ -124,20 +160,24 @@ class Domain(GeometricObject): def __eq__(self, other): """ - Returns true if two sets are equal. + Return True if two domains are equal. """ - symbols = self._xsymbols([self, other]) - islset1 = self._toislset(self.polyhedra, symbols) - islset2 = other._toislset(other.polyhedra, symbols) - equal = bool(libisl.isl_set_is_equal(islset1, islset2)) - libisl.isl_set_free(islset1) - libisl.isl_set_free(islset2) - return equal + if isinstance(other, Domain): + symbols = self._xsymbols([self, other]) + islset1 = self._toislset(self.polyhedra, symbols) + islset2 = other._toislset(other.polyhedra, symbols) + equal = bool(libisl.isl_set_is_equal(islset1, islset2)) + libisl.isl_set_free(islset1) + libisl.isl_set_free(islset2) + return equal + return NotImplemented def isdisjoint(self, other): """ - Return True if two sets have a null intersection. + Return True if two domains have a null intersection. """ + if not isinstance(other, Domain): + raise TypeError('other must be a Domain instance') symbols = self._xsymbols([self, other]) islset1 = self._toislset(self.polyhedra, symbols) islset2 = self._toislset(other.polyhedra, symbols) @@ -148,60 +188,84 @@ class Domain(GeometricObject): def issubset(self, other): """ - Report whether another set contains this set. + Report whether another domain contains the domain. """ - symbols = self._xsymbols([self, other]) - islset1 = self._toislset(self.polyhedra, symbols) - islset2 = self._toislset(other.polyhedra, symbols) - equal = bool(libisl.isl_set_is_subset(islset1, islset2)) - libisl.isl_set_free(islset1) - libisl.isl_set_free(islset2) - return equal + return self < other def __le__(self, other): - """ - Returns true if this set is less than or equal to another set. - """ - return self.issubset(other) + if isinstance(other, Domain): + symbols = self._xsymbols([self, other]) + islset1 = self._toislset(self.polyhedra, symbols) + islset2 = self._toislset(other.polyhedra, symbols) + equal = bool(libisl.isl_set_is_subset(islset1, islset2)) + libisl.isl_set_free(islset1) + libisl.isl_set_free(islset2) + return equal + return NotImplemented + __le__.__doc__ = issubset.__doc__ def __lt__(self, other): """ - Returns true if this set is less than another set. + Report whether another domain is contained within the domain. """ - symbols = self._xsymbols([self, other]) - islset1 = self._toislset(self.polyhedra, symbols) - islset2 = self._toislset(other.polyhedra, symbols) - equal = bool(libisl.isl_set_is_strict_subset(islset1, islset2)) - libisl.isl_set_free(islset1) - libisl.isl_set_free(islset2) - return equal + if isinstance(other, Domain): + symbols = self._xsymbols([self, other]) + islset1 = self._toislset(self.polyhedra, symbols) + islset2 = self._toislset(other.polyhedra, symbols) + equal = bool(libisl.isl_set_is_strict_subset(islset1, islset2)) + libisl.isl_set_free(islset1) + libisl.isl_set_free(islset2) + return equal + return NotImplemented def complement(self): """ - Returns the complement of this set. + Return the complementary domain of the domain. """ islset = self._toislset(self.polyhedra, self.symbols) islset = libisl.isl_set_complement(islset) return self._fromislset(islset, self.symbols) def __invert__(self): + return self.complement() + __invert__.__doc__ = complement.__doc__ + + def make_disjoint(self): """ - Returns the complement of this set. + Return an equivalent domain, whose polyhedra are disjoint. """ - return self.complement() + islset = self._toislset(self.polyhedra, self.symbols) + islset = libisl.isl_set_make_disjoint(islset) + return self._fromislset(islset, self.symbols) - def simplify(self): + def coalesce(self): """ - Returns a set without redundant constraints. + Simplify the representation of the domain by trying to combine pairs of + polyhedra into a single polyhedron, and return the resulting domain. """ islset = self._toislset(self.polyhedra, self.symbols) - islset = libisl.isl_set_remove_redundancies(islset) + islset = libisl.isl_set_coalesce(islset) return self._fromislset(islset, self.symbols) - def aspolyhedron(self): + def detect_equalities(self): + """ + Simplify the representation of the domain by detecting implicit + equalities, and return the resulting domain. """ - Returns polyhedral hull of set. + islset = self._toislset(self.polyhedra, self.symbols) + islset = libisl.isl_set_detect_equalities(islset) + return self._fromislset(islset, self.symbols) + + def remove_redundancies(self): + """ + Remove redundant constraints in the domain, and return the resulting + domain. """ + islset = self._toislset(self.polyhedra, self.symbols) + islset = libisl.isl_set_remove_redundancies(islset) + return self._fromislset(islset, self.symbols) + + def aspolyhedron(self): from .polyhedra import Polyhedron islset = self._toislset(self.polyhedra, self.symbols) islbset = libisl.isl_set_polyhedral_hull(islset) @@ -210,26 +274,33 @@ class Domain(GeometricObject): def asdomain(self): return self - def project(self, dims): + def project(self, symbols): """ - Return new set with given dimensions removed. + Project out the sequence of symbols given in arguments, and return the + resulting domain. """ + symbols = list(symbols) + for symbol in symbols: + if not isinstance(symbol, Symbol): + raise TypeError('symbols must be Symbol instances') islset = self._toislset(self.polyhedra, self.symbols) n = 0 for index, symbol in reversed(list(enumerate(self.symbols))): - if symbol in dims: + if symbol in symbols: n += 1 elif n > 0: - islset = libisl.isl_set_project_out(islset, libisl.isl_dim_set, index + 1, n) + islset = libisl.isl_set_project_out(islset, + libisl.isl_dim_set, index + 1, n) n = 0 if n > 0: islset = libisl.isl_set_project_out(islset, libisl.isl_dim_set, 0, n) - dims = [symbol for symbol in self.symbols if symbol not in dims] - return Domain._fromislset(islset, dims) + symbols = [symbol for symbol in self.symbols if symbol not in symbols] + return Domain._fromislset(islset, symbols) def sample(self): """ - Returns a single subset of the input. + Return a sample of the domain, as an integer instance of Point. If the + domain is empty, a ValueError exception is raised. """ islset = self._toislset(self.polyhedra, self.symbols) islpoint = libisl.isl_set_sample_point(islset) @@ -247,67 +318,67 @@ class Domain(GeometricObject): def intersection(self, *others): """ - Return the intersection of two sets as a new set. + Return the intersection of two or more domains as a new domain. As an + alternative, function And() can be used. """ - if len(others) == 0: - return self - symbols = self._xsymbols((self,) + others) - islset1 = self._toislset(self.polyhedra, symbols) + result = self for other in others: - islset2 = other._toislset(other.polyhedra, symbols) - islset1 = libisl.isl_set_intersect(islset1, islset2) - return self._fromislset(islset1, symbols) + result &= other + return result def __and__(self, other): - """ - Return the intersection of two sets as a new set. - """ - return self.intersection(other) + if isinstance(other, Domain): + symbols = self._xsymbols([self, other]) + islset1 = self._toislset(self.polyhedra, symbols) + islset2 = other._toislset(other.polyhedra, symbols) + islset = libisl.isl_set_intersect(islset1, islset2) + return self._fromislset(islset, symbols) + return NotImplemented + __and__.__doc__ = intersection.__doc__ def union(self, *others): """ - Return the union of sets as a new set. + Return the union of two or more domains as a new domain. As an + alternative, function Or() can be used. """ - if len(others) == 0: - return self - symbols = self._xsymbols((self,) + others) - islset1 = self._toislset(self.polyhedra, symbols) + result = self for other in others: - islset2 = other._toislset(other.polyhedra, symbols) - islset1 = libisl.isl_set_union(islset1, islset2) - return self._fromislset(islset1, symbols) + result |= other + return result def __or__(self, other): - """ - Return a new set with elements from both sets. - """ - return self.union(other) + if isinstance(other, Domain): + symbols = self._xsymbols([self, other]) + islset1 = self._toislset(self.polyhedra, symbols) + islset2 = other._toislset(other.polyhedra, symbols) + islset = libisl.isl_set_union(islset1, islset2) + return self._fromislset(islset, symbols) + return NotImplemented + __or__.__doc__ = union.__doc__ def __add__(self, other): - """ - Return new set containing all elements in both sets. - """ - return self.union(other) + return self | other + __add__.__doc__ = union.__doc__ def difference(self, other): """ - Return the difference of two sets as a new set. + Return the difference of two domains as a new domain. """ - symbols = self._xsymbols([self, other]) - islset1 = self._toislset(self.polyhedra, symbols) - islset2 = other._toislset(other.polyhedra, symbols) - islset = libisl.isl_set_subtract(islset1, islset2) - return self._fromislset(islset, symbols) + return self - other def __sub__(self, other): - """ - Return the difference of two sets as a new set. - """ - return self.difference(other) + if isinstance(other, Domain): + symbols = self._xsymbols([self, other]) + islset1 = self._toislset(self.polyhedra, symbols) + islset2 = other._toislset(other.polyhedra, symbols) + islset = libisl.isl_set_subtract(islset1, islset2) + return self._fromislset(islset, symbols) + return NotImplemented + __sub__.__doc__ = difference.__doc__ def lexmin(self): """ - Return a new set containing the lexicographic minimum of the elements in the set. + Return the lexicographic minimum of the elements in the domain. """ islset = self._toislset(self.polyhedra, self.symbols) islset = libisl.isl_set_lexmin(islset) @@ -315,52 +386,35 @@ class Domain(GeometricObject): def lexmax(self): """ - Return a new set containing the lexicographic maximum of the elements in the set. + Return the lexicographic maximum of the elements in the domain. """ islset = self._toislset(self.polyhedra, self.symbols) islset = libisl.isl_set_lexmax(islset) return self._fromislset(islset, self.symbols) - - def involves_vars(self, vars): - """ - Returns true if a set depends on given dimensions. - """ - islset = self._toislset(self.polyhedra, self.symbols) - dims = sorted(vars) - symbols = sorted(list(self.symbols)) - n = 0 - if len(dims)>0: - for dim in dims: - if dim in symbols: - first = symbols.index(dims[0]) - n +=1 - else: - first = 0 - else: - return False - value = bool(libisl.isl_set_involves_dims(islset, libisl.isl_dim_set, first, n)) - libisl.isl_set_free(islset) - return value - - _RE_COORDINATE = re.compile(r'\((?P\-?\d+)\)(/(?P\d+))?') + if islhelper.isl_version >= '0.13': + _RE_COORDINATE = re.compile(r'\((?P\-?\d+)\)(/(?P\d+))?') + else: + _RE_COORDINATE = None def vertices(self): """ - Return a list of vertices for this Polygon. + Return the vertices of the domain, as a list of rational instances of + Point. """ from .polyhedra import Polyhedron if not self.isbounded(): raise ValueError('domain must be bounded') - islbset = self._toislbasicset(self.equalities, self.inequalities, self.symbols) + islbset = self._toislbasicset(self.equalities, self.inequalities, + self.symbols) vertices = libisl.isl_basic_set_compute_vertices(islbset); vertices = islhelper.isl_vertices_vertices(vertices) points = [] for vertex in vertices: - expr = libisl.isl_vertex_get_expr(vertex) + expression = libisl.isl_vertex_get_expr(vertex) coordinates = [] - if islhelper.isl_version < '0.13': - constraints = islhelper.isl_basic_set_constraints(expr) + if self._RE_COORDINATE is None: + constraints = islhelper.isl_basic_set_constraints(expression) for constraint in constraints: constant = libisl.isl_constraint_get_constant_val(constraint) constant = islhelper.isl_val_to_int(constant) @@ -372,7 +426,7 @@ class Domain(GeometricObject): coordinate = -Fraction(constant, coefficient) coordinates.append((symbol, coordinate)) else: - string = islhelper.isl_multi_aff_to_str(expr) + string = islhelper.isl_multi_aff_to_str(expression) matches = self._RE_COORDINATE.finditer(string) for symbol, match in zip(self.symbols, matches): numerator = int(match.group('num')) @@ -385,7 +439,9 @@ class Domain(GeometricObject): def points(self): """ - Returns the points contained in the set. + Return the integer points of a bounded domain, as a list of integer + instances of Point. If the domain is not bounded, a ValueError exception + is raised. """ if not self.isbounded(): raise ValueError('domain must be bounded') @@ -403,6 +459,15 @@ class Domain(GeometricObject): points.append(Point(coordinates)) return points + def __contains__(self, point): + """ + Return True if the point if contained within the domain. + """ + for polyhedron in self.polyhedra: + if point in polyhedron: + return True + return False + @classmethod def _polygon_inner_point(cls, points): symbols = points[0].symbols @@ -456,7 +521,9 @@ class Domain(GeometricObject): def faces(self): """ - Returns the vertices of the faces of a polyhedra. + Return the list of faces of a bounded domain. Each face is represented + by a list of vertices, in the form of rational instances of Point. If + the domain is not bounded, a ValueError exception is raised. """ faces = [] for polyhedron in self.polyhedra: @@ -518,10 +585,13 @@ class Domain(GeometricObject): axes.set_zlim(zmin, zmax) return axes - def plot(self, plot=None, **kwargs): """ - Display plot of this set. + Plot a 2D or 3D domain using matplotlib. Draw it to the current plot + object if present, otherwise create a new one. options are keyword + arguments passed to the matplotlib drawing functions, they can be used + to set the drawing color for example. Raise ValueError is the domain is + not 2D or 3D. """ if not self.isbounded(): raise ValueError('domain must be bounded') @@ -530,18 +600,14 @@ class Domain(GeometricObject): elif self.dimension == 3: return self._plot_3d(plot=plot, **kwargs) else: - raise ValueError('polyhedron must be 2 or 3-dimensional') - - def __contains__(self, point): - for polyhedron in self.polyhedra: - if point in polyhedron: - return True - return False + raise ValueError('domain must be 2 or 3-dimensional') def subs(self, symbol, expression=None): """ - Subsitute the given value into an expression and return the resulting - expression. + Substitute the given symbol by an expression in the domain constraints. + To perform multiple substitutions at once, pass a sequence or a + dictionary of (old, new) pairs to subs. The syntax of this function is + similar to LinExpr.subs(). """ polyhedra = [polyhedron.subs(symbol, expression) for polyhedron in self.polyhedra] @@ -603,10 +669,10 @@ class Domain(GeometricObject): elif isinstance(node, ast.Compare): equalities = [] inequalities = [] - left = Expression._fromast(node.left) + left = LinExpr._fromast(node.left) for i in range(len(node.ops)): op = node.ops[i] - right = Expression._fromast(node.comparators[i]) + right = LinExpr._fromast(node.comparators[i]) if isinstance(op, ast.Lt): inequalities.append(right - left - 1) elif isinstance(op, ast.LtE): @@ -629,22 +695,26 @@ class Domain(GeometricObject): _RE_AND = re.compile(r'\band\b|,|&&|/\\|∧|∩') _RE_OR = re.compile(r'\bor\b|;|\|\||\\/|∨|∪') _RE_NOT = re.compile(r'\bnot\b|!|¬') - _RE_NUM_VAR = Expression._RE_NUM_VAR + _RE_NUM_VAR = LinExpr._RE_NUM_VAR _RE_OPERATORS = re.compile(r'(&|\||~)') @classmethod def fromstring(cls, string): - # remove curly brackets + """ + Create a domain from a string. Raise SyntaxError if the string is not + properly formatted. + """ + # Remove curly brackets. string = cls._RE_BRACES.sub(r'', string) - # replace '=' by '==' + # Replace '=' by '=='. string = cls._RE_EQ.sub(r'\1==\2', string) - # replace 'and', 'or', 'not' + # Replace 'and', 'or', 'not'. string = cls._RE_AND.sub(r' & ', string) string = cls._RE_OR.sub(r' | ', string) string = cls._RE_NOT.sub(r' ~', string) - # add implicit multiplication operators, e.g. '5x' -> '5*x' + # Add implicit multiplication operators, e.g. '5x' -> '5*x'. string = cls._RE_NUM_VAR.sub(r'\1*\2', string) - # add parentheses to force precedence + # Add parentheses to force precedence. tokens = cls._RE_OPERATORS.split(string) for i, token in enumerate(tokens): if i % 2 == 0: @@ -659,14 +729,11 @@ class Domain(GeometricObject): strings = [repr(polyhedron) for polyhedron in self.polyhedra] return 'Or({})'.format(', '.join(strings)) - def _repr_latex_(self): - strings = [] - for polyhedron in self.polyhedra: - strings.append('({})'.format(polyhedron._repr_latex_().strip('$'))) - return '${}$'.format(' \\vee '.join(strings)) - @classmethod - def fromsympy(cls, expr): + def fromsympy(cls, expression): + """ + Create a domain from a SymPy expression. + """ import sympy from .polyhedra import Lt, Le, Eq, Ne, Ge, Gt funcmap = { @@ -675,14 +742,17 @@ class Domain(GeometricObject): sympy.Eq: Eq, sympy.Ne: Ne, sympy.Ge: Ge, sympy.Gt: Gt, } - if expr.func in funcmap: - args = [Domain.fromsympy(arg) for arg in expr.args] - return funcmap[expr.func](*args) - elif isinstance(expr, sympy.Expr): - return Expression.fromsympy(expr) - raise ValueError('non-domain expression: {!r}'.format(expr)) + if expression.func in funcmap: + args = [Domain.fromsympy(arg) for arg in expression.args] + return funcmap[expression.func](*args) + elif isinstance(expression, sympy.Expr): + return LinExpr.fromsympy(expression) + raise ValueError('non-domain expression: {!r}'.format(expression)) def tosympy(self): + """ + Convert the domain to a SymPy expression. + """ import sympy polyhedra = [polyhedron.tosympy() for polyhedron in polyhedra] return sympy.Or(*polyhedra) @@ -690,7 +760,7 @@ class Domain(GeometricObject): def And(*domains): """ - Return the intersection of two sets as a new set. + Create the intersection domain of the domains given in arguments. """ if len(domains) == 0: from .polyhedra import Universe @@ -700,7 +770,7 @@ def And(*domains): def Or(*domains): """ - Return the union of sets as a new set. + Create the union domain of the domains given in arguments. """ if len(domains) == 0: from .polyhedra import Empty @@ -710,6 +780,6 @@ def Or(*domains): def Not(domain): """ - Returns the complement of this set. + Create the complementary domain of the domain given in argument. """ return ~domain