X-Git-Url: https://scm.cri.ensmp.fr/git/linpy.git/blobdiff_plain/7b93cea1daf2889e9ee10ca9c22a1b5124404937..5a5fd1db359b190c6207301eb08705a34367968a:/linpy/domains.py diff --git a/linpy/domains.py b/linpy/domains.py index 21e78db..83f6c2c 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, Rational 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 and complementary operations. + + A domain with a unique polyhedron is automatically subclassed as a + Polyhedron instance. + """ __slots__ = ( '_polyhedra', @@ -44,6 +51,29 @@ class Domain(GeometricObject): ) def __new__(cls, *polyhedra): + """ + Return a domain from a sequence of polyhedra. + + >>> square = Polyhedron('0 <= x <= 2, 0 <= y <= 2') + >>> square2 = Polyhedron('2 <= x <= 4, 2 <= y <= 4') + >>> dom = Domain([square, square2]) + + It is also possible to build domains from polyhedra using arithmetic + operators Domain.__and__(), Domain.__or__() or functions And() and Or(), + using one of the following instructions: + + >>> square = Polyhedron('0 <= x <= 2, 0 <= y <= 2') + >>> square2 = Polyhedron('2 <= x <= 4, 2 <= y <= 4') + >>> dom = square | square2 + >>> dom = Or(square, square2) + + Alternatively, a domain can be built from a string: + + >>> dom = Domain('0 <= x <= 2, 0 <= y <= 2; 2 <= x <= 4, 2 <= y <= 4') + + 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 +104,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 +134,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 +150,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 +159,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 +187,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(mainctx, 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): """ - Returns polyhedral hull of set. + Simplify the representation of the domain by detecting implicit + equalities, and return the resulting domain. """ + 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 +273,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 +317,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,51 +385,34 @@ 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) coordinates = [] - if islhelper.isl_version < '0.13': + if self._RE_COORDINATE is None: constraints = islhelper.isl_basic_set_constraints(expr) for constraint in constraints: constant = libisl.isl_constraint_get_constant_val(constraint) @@ -385,7 +438,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 +458,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 +520,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 +584,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 +599,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 +668,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,11 +694,15 @@ 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): + """ + 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 '==' @@ -667,6 +736,9 @@ class Domain(GeometricObject): @classmethod def fromsympy(cls, expr): + """ + Create a domain from a sympy expression. + """ import sympy from .polyhedra import Lt, Le, Eq, Ne, Ge, Gt funcmap = { @@ -679,10 +751,13 @@ class Domain(GeometricObject): args = [Domain.fromsympy(arg) for arg in expr.args] return funcmap[expr.func](*args) elif isinstance(expr, sympy.Expr): - return Expression.fromsympy(expr) + return LinExpr.fromsympy(expr) raise ValueError('non-domain expression: {!r}'.format(expr)) def tosympy(self): + """ + Convert the domain to a sympy expression. + """ import sympy polyhedra = [polyhedron.tosympy() for polyhedron in polyhedra] return sympy.Or(*polyhedra) @@ -690,26 +765,29 @@ 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 return Universe else: return domains[0].intersection(*domains[1:]) +And.__doc__ = Domain.intersection.__doc__ 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 return Empty else: return domains[0].union(*domains[1:]) +Or.__doc__ = Domain.union.__doc__ def Not(domain): """ - Returns the complement of this set. + Create the complementary domain of the domain given in argument. """ return ~domain +Not.__doc__ = Domain.complement.__doc__