X-Git-Url: https://scm.cri.ensmp.fr/git/linpy.git/blobdiff_plain/3ff02ec4a397d81d1019080a78198c454861c1e2..9e103619134c8839d143951c5b814bab500fcdcf:/pypol/linear.py diff --git a/pypol/linear.py b/pypol/linear.py index 331d3af..e6d442e 100644 --- a/pypol/linear.py +++ b/pypol/linear.py @@ -1,19 +1,19 @@ -import ctypes, ctypes.util +import ast import functools import numbers +import re from fractions import Fraction, gcd -from . import isl, islhelper -from .isl import libisl, Context +from pypol import isl +from pypol.isl import libisl __all__ = [ - 'Expression', - 'constant', 'symbol', 'symbols', + 'Expression', 'Constant', 'Symbol', 'symbols', 'eq', 'le', 'lt', 'ge', 'gt', 'Polyhedron', - 'empty', 'universe' + 'Empty', 'Universe' ] @@ -23,7 +23,7 @@ def _polymorphic_method(func): if isinstance(b, Expression): return func(a, b) if isinstance(b, numbers.Rational): - b = constant(b) + b = Constant(b) return func(a, b) return NotImplemented return wrapper @@ -34,7 +34,7 @@ def _polymorphic_operator(func): @functools.wraps(func) def wrapper(a, b): if isinstance(a, numbers.Rational): - a = constant(a) + a = Constant(a) return func(a, b) elif isinstance(a, Expression): return func(a, b) @@ -42,48 +42,103 @@ def _polymorphic_operator(func): return wrapper +_main_ctx = isl.Context() + + class Expression: """ This class implements linear expressions. """ + __slots__ = ( + '_coefficients', + '_constant', + '_symbols', + '_dimension', + ) + def __new__(cls, coefficients=None, constant=0): if isinstance(coefficients, str): if constant: raise TypeError('too many arguments') return cls.fromstring(coefficients) - self = super().__new__(cls) - self._coefficients = {} if isinstance(coefficients, dict): coefficients = coefficients.items() - if coefficients is not None: - for symbol, coefficient in coefficients: - if isinstance(symbol, Expression) and symbol.issymbol(): - symbol = str(symbol) - elif not isinstance(symbol, str): - raise TypeError('symbols must be strings') - if not isinstance(coefficient, numbers.Rational): - raise TypeError('coefficients must be rational numbers') - if coefficient != 0: - self._coefficients[symbol] = coefficient + if coefficients is None: + return Constant(constant) + coefficients = [(symbol, coefficient) + for symbol, coefficient in coefficients if coefficient != 0] + if len(coefficients) == 0: + return Constant(constant) + elif len(coefficients) == 1 and constant == 0: + symbol, coefficient = coefficients[0] + if coefficient == 1: + return Symbol(symbol) + self = object().__new__(cls) + self._coefficients = {} + for symbol, coefficient in coefficients: + if isinstance(symbol, Symbol): + symbol = symbol.name + elif not isinstance(symbol, str): + raise TypeError('symbols must be strings or Symbol instances') + if isinstance(coefficient, Constant): + coefficient = coefficient.constant + if not isinstance(coefficient, numbers.Rational): + raise TypeError('coefficients must be rational numbers or Constant instances') + self._coefficients[symbol] = coefficient + if isinstance(constant, Constant): + constant = constant.constant if not isinstance(constant, numbers.Rational): - raise TypeError('constant must be a rational number') + raise TypeError('constant must be a rational number or a Constant instance') self._constant = constant + self._symbols = tuple(sorted(self._coefficients)) + self._dimension = len(self._symbols) return self + @classmethod + def _fromast(cls, node): + if isinstance(node, ast.Module) and len(node.body) == 1: + return cls._fromast(node.body[0]) + elif isinstance(node, ast.Expr): + return cls._fromast(node.value) + elif isinstance(node, ast.Name): + return Symbol(node.id) + elif isinstance(node, ast.Num): + return Constant(node.n) + elif isinstance(node, ast.UnaryOp) and isinstance(node.op, ast.USub): + return -cls._fromast(node.operand) + elif isinstance(node, ast.BinOp): + left = cls._fromast(node.left) + right = cls._fromast(node.right) + if isinstance(node.op, ast.Add): + return left + right + elif isinstance(node.op, ast.Sub): + return left - right + elif isinstance(node.op, ast.Mult): + return left * right + elif isinstance(node.op, ast.Div): + return left / right + raise SyntaxError('invalid syntax') + @classmethod + def fromstring(cls, string): + string = re.sub(r'(\d+|\))\s*([^\W\d_]\w*|\()', r'\1*\2', string) + tree = ast.parse(string, 'eval') + return cls._fromast(tree) + + @property def symbols(self): - yield from sorted(self._coefficients) + return self._symbols @property def dimension(self): - return len(list(self.symbols())) + return self._dimension def coefficient(self, symbol): - if isinstance(symbol, Expression) and symbol.issymbol(): + if isinstance(symbol, Symbol): symbol = str(symbol) elif not isinstance(symbol, str): - raise TypeError('symbol must be a string') + raise TypeError('symbol must be a string or a Symbol instance') try: return self._coefficients[symbol] except KeyError: @@ -91,9 +146,8 @@ class Expression: __getitem__ = coefficient - @property def coefficients(self): - for symbol in self.symbols(): + for symbol in self.symbols: yield symbol, self.coefficient(symbol) @property @@ -101,29 +155,18 @@ class Expression: return self._constant def isconstant(self): - return len(self._coefficients) == 0 + return False def values(self): - for symbol in self.symbols(): + for symbol in self.symbols: yield self.coefficient(symbol) yield self.constant - def values_int(self): - for symbol in self.symbols(): - return self.coefficient(symbol) - return int(self.constant) - - def symbol(self): - if not self.issymbol(): - raise ValueError('not a symbol: {}'.format(self)) - for symbol in self.symbols(): - return symbol - def issymbol(self): - return len(self._coefficients) == 1 and self._constant == 0 + return False def __bool__(self): - return (not self.isconstant()) or bool(self.constant) + return True def __pos__(self): return self @@ -133,8 +176,8 @@ class Expression: @_polymorphic_method def __add__(self, other): - coefficients = dict(self.coefficients) - for symbol, coefficient in other.coefficients: + coefficients = dict(self.coefficients()) + for symbol, coefficient in other.coefficients(): if symbol in coefficients: coefficients[symbol] += coefficient else: @@ -146,8 +189,8 @@ class Expression: @_polymorphic_method def __sub__(self, other): - coefficients = dict(self.coefficients) - for symbol, coefficient in other.coefficients: + coefficients = dict(self.coefficients()) + for symbol, coefficient in other.coefficients(): if symbol in coefficients: coefficients[symbol] -= coefficient else: @@ -161,7 +204,7 @@ class Expression: @_polymorphic_method def __mul__(self, other): if other.isconstant(): - coefficients = dict(self.coefficients) + coefficients = dict(self.coefficients()) for symbol in coefficients: coefficients[symbol] *= other.constant constant = self.constant * other.constant @@ -199,10 +242,9 @@ class Expression: def __str__(self): string = '' - symbols = sorted(self.symbols()) i = 0 - for symbol in symbols: - coefficient = self[symbol] + for symbol in self.symbols: + coefficient = self.coefficient(symbol) if coefficient == 1: if i == 0: string += symbol @@ -244,17 +286,13 @@ class Expression: def __repr__(self): string = '{}({{'.format(self.__class__.__name__) - for i, (symbol, coefficient) in enumerate(self.coefficients): + for i, (symbol, coefficient) in enumerate(self.coefficients()): if i != 0: string += ', ' string += '{!r}: {!r}'.format(symbol, coefficient) string += '}}, {!r})'.format(self.constant) return string - @classmethod - def fromstring(cls, string): - raise NotImplementedError - @_polymorphic_method def __eq__(self, other): # "normal" equality @@ -264,70 +302,119 @@ class Expression: self.constant == other.constant def __hash__(self): - return hash((self._coefficients, self._constant)) + return hash((tuple(sorted(self._coefficients.items())), self._constant)) - def _canonify(self): + def _toint(self): lcm = functools.reduce(lambda a, b: a*b // gcd(a, b), [value.denominator for value in self.values()]) return self * lcm @_polymorphic_method def _eq(self, other): - return Polyhedron(equalities=[(self - other)._canonify()]) + return Polyhedron(equalities=[(self - other)._toint()]) @_polymorphic_method def __le__(self, other): - return Polyhedron(inequalities=[(other - self)._canonify()]) + return Polyhedron(inequalities=[(other - self)._toint()]) @_polymorphic_method def __lt__(self, other): - return Polyhedron(inequalities=[(other - self)._canonify() - 1]) + return Polyhedron(inequalities=[(other - self)._toint() - 1]) @_polymorphic_method def __ge__(self, other): - return Polyhedron(inequalities=[(self - other)._canonify()]) + return Polyhedron(inequalities=[(self - other)._toint()]) @_polymorphic_method def __gt__(self, other): - return Polyhedron(inequalities=[(self - other)._canonify() - 1]) + return Polyhedron(inequalities=[(self - other)._toint() - 1]) + + +class Constant(Expression): + + def __new__(cls, numerator=0, denominator=None): + self = object().__new__(cls) + if denominator is None: + if isinstance(numerator, numbers.Rational): + self._constant = numerator + elif isinstance(numerator, Constant): + self._constant = numerator.constant + else: + raise TypeError('constant must be a rational number or a Constant instance') + else: + self._constant = Fraction(numerator, denominator) + self._coefficients = {} + self._symbols = () + self._dimension = 0 + return self + + def isconstant(self): + return True + def __bool__(self): + return bool(self.constant) -def constant(numerator=0, denominator=None): - if denominator is None and isinstance(numerator, numbers.Rational): - return Expression(constant=numerator) - else: - return Expression(constant=Fraction(numerator, denominator)) + def __repr__(self): + if self.constant.denominator == 1: + return '{}({!r})'.format(self.__class__.__name__, self.constant) + else: + return '{}({!r}, {!r})'.format(self.__class__.__name__, + self.constant.numerator, self.constant.denominator) + +class Symbol(Expression): + + __slots__ = Expression.__slots__ + ( + '_name', + ) + + def __new__(cls, name): + if isinstance(name, Symbol): + name = name.name + elif not isinstance(name, str): + raise TypeError('name must be a string or a Symbol instance') + self = object().__new__(cls) + self._coefficients = {name: 1} + self._constant = 0 + self._symbols = tuple(name) + self._name = name + self._dimension = 1 + return self -def symbol(name): - if not isinstance(name, str): - raise TypeError('name must be a string') - return Expression(coefficients={name: 1}) + @property + def name(self): + return self._name + + def issymbol(self): + return True + + def __repr__(self): + return '{}({!r})'.format(self.__class__.__name__, self._name) def symbols(names): if isinstance(names, str): names = names.replace(',', ' ').split() - return (symbol(name) for name in names) + return (Symbol(name) for name in names) @_polymorphic_operator def eq(a, b): - return a._eq(b) + return a.__eq__(b) @_polymorphic_operator def le(a, b): - return a <= b + return a.__le__(b) @_polymorphic_operator def lt(a, b): - return a < b + return a.__lt__(b) @_polymorphic_operator def ge(a, b): - return a >= b + return a.__ge__(b) @_polymorphic_operator def gt(a, b): - return a > b + return a.__gt__(b) class Polyhedron: @@ -335,6 +422,13 @@ class Polyhedron: This class implements polyhedrons. """ + __slots__ = ( + '_equalities', + '_inequalities', + '_constraints', + '_symbols', + ) + def __new__(cls, equalities=None, inequalities=None): if isinstance(equalities, str): if inequalities is not None: @@ -349,6 +443,7 @@ class Polyhedron: raise TypeError('non-integer constraint: ' '{} == 0'.format(constraint)) self._equalities.append(constraint) + self._equalities = tuple(self._equalities) self._inequalities = [] if inequalities is not None: for constraint in inequalities: @@ -357,92 +452,140 @@ class Polyhedron: raise TypeError('non-integer constraint: ' '{} <= 0'.format(constraint)) self._inequalities.append(constraint) - self._bset = self._to_isl() - #print(self._bset) - #put this here just to test from isl method - #from_isl = self.from_isl(self._bset) - #print(from_isl) - #rint(self) + self._inequalities = tuple(self._inequalities) + self._constraints = self._equalities + self._inequalities + self._symbols = set() + for constraint in self._constraints: + self.symbols.update(constraint.symbols) + self._symbols = tuple(sorted(self._symbols)) return self + @classmethod + def _fromast(cls, node): + if isinstance(node, ast.Module) and len(node.body) == 1: + return cls._fromast(node.body[0]) + elif isinstance(node, ast.Expr): + return cls._fromast(node.value) + elif isinstance(node, ast.BinOp) and isinstance(node.op, ast.BitAnd): + equalities1, inequalities1 = cls._fromast(node.left) + equalities2, inequalities2 = cls._fromast(node.right) + equalities = equalities1 + equalities2 + inequalities = inequalities1 + inequalities2 + return equalities, inequalities + elif isinstance(node, ast.Compare): + equalities = [] + inequalities = [] + left = Expression._fromast(node.left) + for i in range(len(node.ops)): + op = node.ops[i] + right = Expression._fromast(node.comparators[i]) + if isinstance(op, ast.Lt): + inequalities.append(right - left - 1) + elif isinstance(op, ast.LtE): + inequalities.append(right - left) + elif isinstance(op, ast.Eq): + equalities.append(left - right) + elif isinstance(op, ast.GtE): + inequalities.append(left - right) + elif isinstance(op, ast.Gt): + inequalities.append(left - right - 1) + else: + break + left = right + else: + return equalities, inequalities + raise SyntaxError('invalid syntax') + + @classmethod + def fromstring(cls, string): + string = string.strip() + string = re.sub(r'^\{\s*|\s*\}\$', '', string) + string = re.sub(r'([^<=>])=([^<=>])', r'\1==\2', string) + string = re.sub(r'(\d+|\))\s*([^\W\d_]\w*|\()', r'\1*\2', string) + tokens = re.split(r',|;|and|&&|/\\|â§', string, flags=re.I) + tokens = ['({})'.format(token) for token in tokens] + string = ' & '.join(tokens) + tree = ast.parse(string, 'eval') + equalities, inequalities = cls._fromast(tree) + return cls(equalities, inequalities) + @property def equalities(self): - yield from self._equalities + return self._equalities @property def inequalities(self): - yield from self._inequalities + return self._inequalities @property - def constant(self): - return self._constant - - def isconstant(self): - return len(self._coefficients) == 0 - - def isempty(self): - return bool(libisl.isl_basic_set_is_empty(self._bset)) - def constraints(self): - yield from self.equalities - yield from self.inequalities + return self._constraints + @property def symbols(self): - s = set() - for constraint in self.constraints(): - s.update(constraint.symbols()) - return sorted(s) + return self._symbols @property def dimension(self): - return len(self.symbols()) + return len(self.symbols) def __bool__(self): - # return false if the polyhedron is empty, true otherwise - if self._equalities or self._inequalities: - return False - else: - return True + return not self.is_empty() def __contains__(self, value): # is the value in the polyhedron? raise NotImplementedError def __eq__(self, other): - raise NotImplementedError + # works correctly when symbols is not passed + # should be equal if values are the same even if symbols are different + bset = self._toisl() + other = other._toisl() + return bool(libisl.isl_basic_set_plain_is_equal(bset, other)) - def is_empty(self): - return + def isempty(self): + bset = self._toisl() + return bool(libisl.isl_basic_set_is_empty(bset)) def isuniverse(self): - return self == universe + bset = self._toisl() + return bool(libisl.isl_basic_set_is_universe(bset)) def isdisjoint(self, other): # return true if the polyhedron has no elements in common with other - raise NotImplementedError + #symbols = self._symbolunion(other) + bset = self._toisl() + other = other._toisl() + return bool(libisl.isl_set_is_disjoint(bset, other)) def issubset(self, other): - raise NotImplementedError + # check if self(bset) is a subset of other + symbols = self._symbolunion(other) + bset = self._toisl(symbols) + other = other._toisl(symbols) + return bool(libisl.isl_set_is_strict_subset(other, bset)) def __le__(self, other): return self.issubset(other) def __lt__(self, other): - raise NotImplementedError + symbols = self._symbolunion(other) + bset = self._toisl(symbols) + other = other._toisl(symbols) + return bool(libisl.isl_set_is_strict_subset(other, bset)) def issuperset(self, other): # test whether every element in other is in the polyhedron - for value in other: - if value == self.constraints(): - return True - else: - return False raise NotImplementedError def __ge__(self, other): return self.issuperset(other) def __gt__(self, other): + symbols = self._symbolunion(other) + bset = self._toisl(symbols) + other = other._toisl(symbols) + bool(libisl.isl_set_is_strict_subset(other, bset)) raise NotImplementedError def union(self, *others): @@ -468,10 +611,13 @@ class Polyhedron: def __and__(self, other): return self.intersection(other) - def difference(self, *others): - # return a new polyhedron with elements in the polyhedron that are not - # in the others - raise NotImplementedError + def difference(self, other): + # return a new polyhedron with elements in the polyhedron that are not in the other + symbols = self._symbolunion(other) + bset = self._toisl(symbols) + other = other._toisl(symbols) + difference = libisl.isl_set_subtract(bset, other) + return difference def __sub__(self, other): return self.difference(other) @@ -485,83 +631,75 @@ class Polyhedron: return '{{{}}}'.format(', '.join(constraints)) def __repr__(self): - equalities = list(self.equalities) - inequalities = list(self.inequalities) - return '{}(equalities={!r}, inequalities={!r})' \ - ''.format(self.__class__.__name__, equalities, inequalities) - - @classmethod - def fromstring(cls, string): - raise NotImplementedError + if self.isempty(): + return 'Empty' + elif self.isuniverse(): + return 'Universe' + else: + equalities = list(self.equalities) + inequalities = list(self.inequalities) + return '{}(equalities={!r}, inequalities={!r})' \ + ''.format(self.__class__.__name__, equalities, inequalities) def _symbolunion(self, *others): - symbols = set(self.symbols()) + symbols = set(self.symbols) for other in others: - symbols.update(other.symbols()) + symbols.update(other.symbols) return sorted(symbols) - def _to_isl(self, symbols=None): + def _toisl(self, symbols=None): if symbols is None: - symbols = self.symbols() - num_coefficients = len(symbols) - ctx = Context() - space = libisl.isl_space_set_alloc(ctx, 0, num_coefficients) + symbols = self.symbols + dimension = len(symbols) + space = libisl.isl_space_set_alloc(_main_ctx, 0, dimension) bset = libisl.isl_basic_set_universe(libisl.isl_space_copy(space)) ls = libisl.isl_local_space_from_space(space) - ceq = libisl.isl_equality_alloc(libisl.isl_local_space_copy(ls)) - cin = libisl.isl_inequality_alloc(libisl.isl_local_space_copy(ls)) - '''if there are equalities/inequalities, take each constant and coefficient and add as a constraint to the basic set''' - if list(self.equalities): #check if any equalities exist - for eq in self.equalities: - coeff_eq = dict(eq.coefficients) - if eq.constant: - value = eq.constant - ceq = libisl.isl_constraint_set_constant_si(ceq, value) - for eq in coeff_eq: - num = coeff_eq.get(eq) - iden = symbols.index(eq) - ceq = libisl.isl_constraint_set_coefficient_si(ceq, islhelper.isl_dim_set, iden, num) #use 3 for type isl_dim_set + for equality in self.equalities: + ceq = libisl.isl_equality_alloc(libisl.isl_local_space_copy(ls)) + for symbol, coefficient in equality.coefficients(): + val = str(coefficient).encode() + val = libisl.isl_val_read_from_str(_main_ctx, val) + dim = symbols.index(symbol) + ceq = libisl.isl_constraint_set_coefficient_val(ceq, libisl.isl_dim_set, dim, val) + if equality.constant != 0: + val = str(equality.constant).encode() + val = libisl.isl_val_read_from_str(_main_ctx, val) + ceq = libisl.isl_constraint_set_constant_val(ceq, val) bset = libisl.isl_basic_set_add_constraint(bset, ceq) - if list(self.inequalities): #check if any inequalities exist - for ineq in self.inequalities: - coeff_in = dict(ineq.coefficients) - if ineq.constant: - value = ineq.constant - cin = libisl.isl_constraint_set_constant_si(cin, value) - for ineq in coeff_in: - num = coeff_in.get(ineq) - iden = symbols.index(ineq) - cin = libisl.isl_constraint_set_coefficient_si(cin, islhelper.isl_dim_set, iden, num) #use 3 for type isl_dim_set + for inequality in self.inequalities: + cin = libisl.isl_inequality_alloc(libisl.isl_local_space_copy(ls)) + for symbol, coefficient in inequality.coefficients(): + val = str(coefficient).encode() + val = libisl.isl_val_read_from_str(_main_ctx, val) + dim = symbols.index(symbol) + cin = libisl.isl_constraint_set_coefficient_val(cin, libisl.isl_dim_set, dim, val) + if inequality.constant != 0: + val = str(inequality.constant).encode() + val = libisl.isl_val_read_from_str(_main_ctx, val) + cin = libisl.isl_constraint_set_constant_val(cin, val) bset = libisl.isl_basic_set_add_constraint(bset, cin) - ip = libisl.isl_printer_to_str(ctx) #create string printer - ip = libisl.isl_printer_print_basic_set(ip, bset) #print basic set to printer - string = libisl.isl_printer_get_str(ip) #get string from printer - string = str(string.decode()) - print(string) + bset = isl.BasicSet(bset) return bset - def from_isl(self, bset): - '''takes basic set in isl form and puts back into python version of polyhedron + @classmethod + def _fromisl(cls, bset, symbols): + raise NotImplementedError + equalities = ... + inequalities = ... + return cls(equalities, inequalities) + '''takes basic set in isl form and puts back into python version of polyhedron isl example code gives isl form as: "{[i] : exists (a : i = 2a and i >= 10 and i <= 42)}") our printer is giving form as: - b'{ [i0] : 1 = 0 }' ''' - #bset = self - if self._equalities: - constraints = libisl.isl_basic_set_equalities_matrix(bset, 3) - elif self._inequalities: - constraints = libisl.isl_basic_set_inequalities_matrix(bset, 3) - print(constraints) - return constraints - -empty = None #eq(0,1) -universe = None #Polyhedron() + { [i0, i1] : 2i1 >= -2 - i0 } ''' -if __name__ == '__main__': - ex1 = Expression(coefficients={'a': 1, 'x': 2}, constant=2) - ex2 = Expression(coefficients={'a': 3 , 'b': 2}, constant=3) - p = Polyhedron(inequalities=[ex1, ex2]) - #p = eq(ex2, 0)# 2a+4 = 0, in fact 6a+3 = 0 - #p.to_isl() +Empty = eq(0,1) + +Universe = Polyhedron() -#universe = Polyhedron() + +if __name__ == '__main__': + p1 = Polyhedron('2a + 2b + 1 == 0') # empty + print(p1._toisl()) + p2 = Polyhedron('3x + 2y + 3 == 0') # not empty + print(p2._toisl())