X-Git-Url: https://svn.cri.ensmp.fr/git/linpy.git/blobdiff_plain/50f0fde4dca3f10e535314b1c9d325cd9283735b..3ffceb83a27c41c65afe87e4370ad5f07ba0a093:/linpy/domains.py diff --git a/linpy/domains.py b/linpy/domains.py index 70e5e4f..519059d 100644 --- a/linpy/domains.py +++ b/linpy/domains.py @@ -15,13 +15,6 @@ # You should have received a copy of the GNU General Public License # along with LinPy. If not, see . -""" -Polyhedral domains - -This module provides classes and functions to deal with polyhedral -domains, i.e. unions of polyhedra. -""" - import ast import functools import re @@ -31,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 @@ -44,7 +37,11 @@ __all__ = [ @functools.total_ordering class Domain(GeometricObject): """ - This class represents polyhedral domains, i.e. unions of polyhedra. + 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__ = ( @@ -55,7 +52,27 @@ class Domain(GeometricObject): def __new__(cls, *polyhedra): """ - Create and return a new domain from a string or a list of 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: @@ -88,35 +105,28 @@ class Domain(GeometricObject): @property def polyhedra(self): """ - The tuple of polyhedra which constitute the domain. + The tuple of polyhedra present in the domain. """ return self._polyhedra @property def symbols(self): """ - The tuple of symbols present in the domain equations. + The tuple of symbols present in the domain equations, sorted according + to Symbol.sortkey(). """ return self._symbols @property def dimension(self): """ - The dimension of the domain, i.e. the number of symbols. + The dimension of the domain, i.e. the number of symbols present in it. """ return self._dimension - def make_disjoint(self): - """ - Return an equivalent domain, whose polyhedra are disjoint. - """ - islset = self._toislset(self.polyhedra, self.symbols) - islset = libisl.isl_set_make_disjoint(mainctx, islset) - return self._fromislset(islset, self.symbols) - def isempty(self): """ - Return True if the domain is empty. + 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)) @@ -131,7 +141,7 @@ class Domain(GeometricObject): def isuniverse(self): """ - Return True if the domain is universal, i.e. with no constraint. + 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)) @@ -149,20 +159,24 @@ class Domain(GeometricObject): def __eq__(self, other): """ - Return True if the two domains 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 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) @@ -175,29 +189,33 @@ class Domain(GeometricObject): """ 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): - 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): """ 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): """ @@ -211,10 +229,18 @@ class Domain(GeometricObject): return self.complement() __invert__.__doc__ = complement.__doc__ + def make_disjoint(self): + """ + Return an equivalent domain, whose polyhedra are disjoint. + """ + islset = self._toislset(self.polyhedra, self.symbols) + islset = libisl.isl_set_make_disjoint(mainctx, islset) + return self._fromislset(islset, self.symbols) + def coalesce(self): """ Simplify the representation of the domain by trying to combine pairs of - polyhedra into a single polyhedron. + polyhedra into a single polyhedron, and return the resulting domain. """ islset = self._toislset(self.polyhedra, self.symbols) islset = libisl.isl_set_coalesce(islset) @@ -223,7 +249,7 @@ class Domain(GeometricObject): def detect_equalities(self): """ Simplify the representation of the domain by detecting implicit - equalities. + equalities, and return the resulting domain. """ islset = self._toislset(self.polyhedra, self.symbols) islset = libisl.isl_set_detect_equalities(islset) @@ -231,16 +257,14 @@ class Domain(GeometricObject): def remove_redundancies(self): """ - Remove redundant constraints in the domain. + 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): - """ - Return the polyhedral hull of the domain. - """ from .polyhedra import Polyhedron islset = self._toislset(self.polyhedra, self.symbols) islbset = libisl.isl_set_polyhedral_hull(islset) @@ -251,8 +275,13 @@ class Domain(GeometricObject): def project(self, symbols): """ - Project out the symbols given in arguments. + 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))): @@ -269,7 +298,8 @@ class Domain(GeometricObject): def sample(self): """ - Return a sample of the domain. + 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) @@ -287,54 +317,62 @@ class Domain(GeometricObject): def intersection(self, *others): """ - Return the intersection of two domains as a new domain. + 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 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 two domains as a new domain. + 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 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 self.union(other) + return self | other __add__.__doc__ = union.__doc__ def difference(self, other): """ 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 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): @@ -353,11 +391,15 @@ class Domain(GeometricObject): islset = libisl.isl_set_lexmax(islset) return self._fromislset(islset, self.symbols) - _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 the vertices of the domain. + Return the vertices of the domain, as a list of rational instances of + Point. """ from .polyhedra import Polyhedron if not self.isbounded(): @@ -370,7 +412,7 @@ class Domain(GeometricObject): 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) @@ -396,7 +438,9 @@ class Domain(GeometricObject): def points(self): """ - Return the points with integer coordinates contained in the domain. + 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') @@ -414,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 @@ -467,7 +520,9 @@ class Domain(GeometricObject): def faces(self): """ - Return the vertices of the domain, grouped by face. + 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: @@ -531,7 +586,11 @@ class Domain(GeometricObject): def plot(self, plot=None, **kwargs): """ - Plot the domain using matplotlib. + 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') @@ -540,21 +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): - """ - Return True if point if contained within the domain. - """ - 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 symbol by expression in equations and return the resulting - domain. + 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] @@ -616,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): @@ -642,25 +694,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): """ - Convert a string into a domain. + Create a domain from a string. Raise SyntaxError if the string is not + properly formatted. """ - # remove curly brackets + # 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: @@ -671,9 +724,6 @@ class Domain(GeometricObject): return cls._fromast(tree) def __repr__(self): - """ - Return repr(self). - """ assert len(self.polyhedra) >= 2 strings = [repr(polyhedron) for polyhedron in self.polyhedra] return 'Or({})'.format(', '.join(strings)) @@ -687,7 +737,7 @@ class Domain(GeometricObject): @classmethod def fromsympy(cls, expr): """ - Convert a SymPy expression into a domain. + Create a domain from a SymPy expression. """ import sympy from .polyhedra import Lt, Le, Eq, Ne, Ge, Gt @@ -701,12 +751,12 @@ 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 a domain into a SymPy expression. + Convert the domain to a SymPy expression. """ import sympy polyhedra = [polyhedron.tosympy() for polyhedron in polyhedra] @@ -714,21 +764,27 @@ class Domain(GeometricObject): def And(*domains): + """ + 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): + """ + 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): + """ + Create the complementary domain of the domain given in argument. + """ return ~domain -Not.__doc__ = Domain.complement.__doc__