pypol/linear.py

   1 import functools
   2 import numbers
   3 import ctypes, ctypes.util
   4 from pypol import isl
   5
   6 from fractions import Fraction, gcd
   7
   8 libisl = ctypes.CDLL(ctypes.util.find_library('isl'))
   9
  10 libisl.isl_printer_get_str.restype = ctypes.c_char_p
  11
  12 __all__ = [
  13     'Expression',
  14     'constant', 'symbol', 'symbols',
  15     'eq', 'le', 'lt', 'ge', 'gt',
  16     'Polyhedron',
  17     'empty', 'universe'
  18 ]
  19
  20
  21 _CONTEXT = isl.Context()
  22
  23 def _polymorphic_method(func):
  24     @functools.wraps(func)
  25     def wrapper(a, b):
  26         if isinstance(b, Expression):
  27             return func(a, b)
  28         if isinstance(b, numbers.Rational):
  29             b = constant(b)
  30             return func(a, b)
  31         return NotImplemented
  32     return wrapper
  33
  34 def _polymorphic_operator(func):
  35     # A polymorphic operator should call a polymorphic method, hence we just
  36     # have to test the left operand.
  37     @functools.wraps(func)
  38     def wrapper(a, b):
  39         if isinstance(a, numbers.Rational):
  40             a = constant(a)
  41             return func(a, b)
  42         elif isinstance(a, Expression):
  43             return func(a, b)
  44         raise TypeError('arguments must be linear expressions')
  45     return wrapper
  46
  47
  48 class Expression:
  49     """
  50     This class implements linear expressions.
  51     """
  52
  53     def __new__(cls, coefficients=None, constant=0):
  54         if isinstance(coefficients, str):
  55             if constant:
  56                 raise TypeError('too many arguments')
  57             return cls.fromstring(coefficients)
  58         self = super().__new__(cls)
  59         self._coefficients = {}
  60         if isinstance(coefficients, dict):
  61             coefficients = coefficients.items()
  62         if coefficients is not None:
  63             for symbol, coefficient in coefficients:
  64                 if isinstance(symbol, Expression) and symbol.issymbol():
  65                     symbol = str(symbol)
  66                 elif not isinstance(symbol, str):
  67                     raise TypeError('symbols must be strings')
  68                 if not isinstance(coefficient, numbers.Rational):
  69                     raise TypeError('coefficients must be rational numbers')
  70                 if coefficient != 0:
  71                     self._coefficients[symbol] = coefficient
  72         if not isinstance(constant, numbers.Rational):
  73             raise TypeError('constant must be a rational number')
  74         self._constant = constant
  75         return self
  76
  77
  78     def symbols(self):
  79         yield from sorted(self._coefficients)
  80
  81     @property
  82     def dimension(self):
  83         return len(list(self.symbols()))
  84
  85     def coefficient(self, symbol):
  86         if isinstance(symbol, Expression) and symbol.issymbol():
  87             symbol = str(symbol)
  88         elif not isinstance(symbol, str):
  89             raise TypeError('symbol must be a string')
  90         try:
  91             return self._coefficients[symbol]
  92         except KeyError:
  93             return 0
  94
  95     __getitem__ = coefficient
  96
  97     def coefficients(self):
  98         for symbol in self.symbols():
  99             yield symbol, self.coefficient(symbol)
 100
 101     @property
 102     def constant(self):
 103         return self._constant
 104
 105     def isconstant(self):
 106         return len(self._coefficients) == 0
 107
 108     def values(self):
 109         for symbol in self.symbols():
 110             yield self.coefficient(symbol)
 111         yield self.constant
 112
 113     def values_int(self):
 114         for symbol in self.symbols():
 115             return self.coefficient(symbol)
 116         return int(self.constant)
 117
 118
 119     def symbol(self):
 120         if not self.issymbol():
 121             raise ValueError('not a symbol: {}'.format(self))
 122         for symbol in self.symbols():
 123             return symbol
 124
 125     def issymbol(self):
 126         return len(self._coefficients) == 1 and self._constant == 0
 127
 128     def __bool__(self):
 129         return (not self.isconstant()) or bool(self.constant)
 130
 131     def __pos__(self):
 132         return self
 133
 134     def __neg__(self):
 135         return self * -1
 136
 137     @_polymorphic_method
 138     def __add__(self, other):
 139         coefficients = dict(self.coefficients())
 140         for symbol, coefficient in other.coefficients():
 141             if symbol in coefficients:
 142                 coefficients[symbol] += coefficient
 143             else:
 144                 coefficients[symbol] = coefficient
 145         constant = self.constant + other.constant
 146         return Expression(coefficients, constant)
 147
 148     __radd__ = __add__
 149
 150     @_polymorphic_method
 151     def __sub__(self, other):
 152         coefficients = dict(self.coefficients())
 153         for symbol, coefficient in other.coefficients():
 154             if symbol in coefficients:
 155                 coefficients[symbol] -= coefficient
 156             else:
 157                 coefficients[symbol] = -coefficient
 158         constant = self.constant - other.constant
 159         return Expression(coefficients, constant)
 160
 161     def __rsub__(self, other):
 162         return -(self - other)
 163
 164     @_polymorphic_method
 165     def __mul__(self, other):
 166         if other.isconstant():
 167             coefficients = dict(self.coefficients())
 168             for symbol in coefficients:
 169                 coefficients[symbol] *= other.constant
 170             constant = self.constant * other.constant
 171             return Expression(coefficients, constant)
 172         if isinstance(other, Expression) and not self.isconstant():
 173             raise ValueError('non-linear expression: '
 174                     '{} * {}'.format(self._parenstr(), other._parenstr()))
 175         return NotImplemented
 176
 177     __rmul__ = __mul__
 178
 179     @_polymorphic_method
 180     def __truediv__(self, other):
 181         if other.isconstant():
 182             coefficients = dict(self.coefficients())
 183             for symbol in coefficients:
 184                 coefficients[symbol] = \
 185                         Fraction(coefficients[symbol], other.constant)
 186             constant = Fraction(self.constant, other.constant)
 187             return Expression(coefficients, constant)
 188         if isinstance(other, Expression):
 189             raise ValueError('non-linear expression: '
 190                 '{} / {}'.format(self._parenstr(), other._parenstr()))
 191         return NotImplemented
 192
 193     def __rtruediv__(self, other):
 194         if isinstance(other, self):
 195             if self.isconstant():
 196                 constant = Fraction(other, self.constant)
 197                 return Expression(constant=constant)
 198             else:
 199                 raise ValueError('non-linear expression: '
 200                         '{} / {}'.format(other._parenstr(), self._parenstr()))
 201         return NotImplemented
 202
 203     def __str__(self):
 204         string = ''
 205         symbols = sorted(self.symbols())
 206         i = 0
 207         for symbol in symbols:
 208             coefficient = self[symbol]
 209             if coefficient == 1:
 210                 if i == 0:
 211                     string += symbol
 212                 else:
 213                     string += ' + {}'.format(symbol)
 214             elif coefficient == -1:
 215                 if i == 0:
 216                     string += '-{}'.format(symbol)
 217                 else:
 218                     string += ' - {}'.format(symbol)
 219             else:
 220                 if i == 0:
 221                     string += '{}*{}'.format(coefficient, symbol)
 222                 elif coefficient > 0:
 223                     string += ' + {}*{}'.format(coefficient, symbol)
 224                 else:
 225                     assert coefficient < 0
 226                     coefficient *= -1
 227                     string += ' - {}*{}'.format(coefficient, symbol)
 228             i += 1
 229         constant = self.constant
 230         if constant != 0 and i == 0:
 231             string += '{}'.format(constant)
 232         elif constant > 0:
 233             string += ' + {}'.format(constant)
 234         elif constant < 0:
 235             constant *= -1
 236             string += ' - {}'.format(constant)
 237         if string == '':
 238             string = '0'
 239         return string
 240
 241     def _parenstr(self, always=False):
 242         string = str(self)
 243         if not always and (self.isconstant() or self.issymbol()):
 244             return string
 245         else:
 246             return '({})'.format(string)
 247
 248     def __repr__(self):
 249         string = '{}({{'.format(self.__class__.__name__)
 250         for i, (symbol, coefficient) in enumerate(self.coefficients()):
 251             if i != 0:
 252                 string += ', '
 253             string += '{!r}: {!r}'.format(symbol, coefficient)
 254         string += '}}, {!r})'.format(self.constant)
 255         return string
 256
 257     @classmethod
 258     def fromstring(cls, string):
 259         raise NotImplementedError
 260
 261     @_polymorphic_method
 262     def __eq__(self, other):
 263         # "normal" equality
 264         # see http://docs.sympy.org/dev/tutorial/gotchas.html#equals-signs
 265         return isinstance(other, Expression) and \
 266                 self._coefficients == other._coefficients and \
 267                 self.constant == other.constant
 268
 269     def __hash__(self):
 270         return hash((self._coefficients, self._constant))
 271
 272     def _canonify(self):
 273         lcm = functools.reduce(lambda a, b: a*b // gcd(a, b),
 274                 [value.denominator for value in self.values()])
 275         return self * lcm
 276
 277     @_polymorphic_method
 278     def _eq(self, other):
 279         return Polyhedron(equalities=[(self - other)._canonify()])
 280
 281     @_polymorphic_method
 282     def __le__(self, other):
 283         return Polyhedron(inequalities=[(self - other)._canonify()])
 284
 285     @_polymorphic_method
 286     def __lt__(self, other):
 287         return Polyhedron(inequalities=[(self - other)._canonify() + 1])
 288
 289     @_polymorphic_method
 290     def __ge__(self, other):
 291         return Polyhedron(inequalities=[(other - self)._canonify()])
 292
 293     @_polymorphic_method
 294     def __gt__(self, other):
 295         return Polyhedron(inequalities=[(other - self)._canonify() + 1])
 296
 297
 298 def constant(numerator=0, denominator=None):
 299     if denominator is None and isinstance(numerator, numbers.Rational):
 300         return Expression(constant=3)
 301     else:
 302         return Expression(constant=Fraction(numerator, denominator))
 303
 304 def symbol(name):
 305     if not isinstance(name, str):
 306         raise TypeError('name must be a string')
 307     return Expression(coefficients={name: 1})
 308
 309 def symbols(names):
 310     if isinstance(names, str):
 311         names = names.replace(',', ' ').split()
 312     return (symbol(name) for name in names)
 313
 314
 315 @_polymorphic_operator
 316 def eq(a, b):
 317     return a._eq(b)
 318
 319 @_polymorphic_operator
 320 def le(a, b):
 321     return a <= b
 322
 323 @_polymorphic_operator
 324 def lt(a, b):
 325     return a < b
 326
 327 @_polymorphic_operator
 328 def ge(a, b):
 329     return a >= b
 330
 331 @_polymorphic_operator
 332 def gt(a, b):
 333     return a > b
 334
 335
 336 class Polyhedron:
 337     """
 338     This class implements polyhedrons.
 339     """
 340
 341     def __new__(cls, equalities=None, inequalities=None):
 342         if isinstance(equalities, str):
 343             if inequalities is not None:
 344                 raise TypeError('too many arguments')
 345             return cls.fromstring(equalities)
 346         self = super().__new__(cls)
 347         self._equalities = []
 348         if equalities is not None:
 349             for constraint in equalities:
 350                 for value in constraint.values():
 351                     if value.denominator != 1:
 352                         raise TypeError('non-integer constraint: '
 353                                 '{} == 0'.format(constraint))
 354                 self._equalities.append(constraint)
 355         self._inequalities = []
 356         if inequalities is not None:
 357             for constraint in inequalities:
 358                 for value in constraint.values():
 359                     if value.denominator != 1:
 360                         raise TypeError('non-integer constraint: '
 361                                 '{} <= 0'.format(constraint))
 362                 self._inequalities.append(constraint)
 363         self._bset = self.to_isl()
 364         return self._bset
 365
 366
 367     @property
 368     def equalities(self):
 369         yield from self._equalities
 370
 371     @property
 372     def inequalities(self):
 373         yield from self._inequalities
 374
 375     @property
 376     def constant(self):
 377         return self._constant
 378
 379     def isconstant(self):
 380         return len(self._coefficients) == 0
 381
 382
 383     def isempty(self):
 384         return bool(libisl.isl_basic_set_is_empty(self._bset))
 385
 386     def constraints(self):
 387         yield from self.equalities
 388         yield from self.inequalities
 389
 390
 391     def symbols(self):
 392         s = set()
 393         for constraint in self.constraints():
 394             s.update(constraint.symbols)
 395             yield from sorted(s)
 396
 397     def symbol_count(self):
 398         s = []
 399         for constraint in self.constraints():
 400             s.append(constraint.symbols)
 401         return s
 402
 403     @property
 404     def dimension(self):
 405         return len(self.symbols())
 406
 407     def __bool__(self):
 408         # return false if the polyhedron is empty, true otherwise
 409         if self._equalities or self._inequalities:
 410             return False
 411         else:
 412             return True
 413
 414
 415     def __contains__(self, value):
 416         # is the value in the polyhedron?
 417         raise NotImplementedError
 418
 419     def __eq__(self, other):
 420         raise NotImplementedError
 421
 422     def is_empty(self):
 423         return
 424
 425     def isuniverse(self):
 426         return self == universe
 427
 428     def isdisjoint(self, other):
 429         # return true if the polyhedron has no elements in common with other
 430         raise NotImplementedError
 431
 432     def issubset(self, other):
 433         raise NotImplementedError
 434
 435     def __le__(self, other):
 436         return self.issubset(other)
 437
 438     def __lt__(self, other):
 439         raise NotImplementedError
 440
 441     def issuperset(self, other):
 442         # test whether every element in other is in the polyhedron
 443         for value in other:
 444             if value == self.constraints():
 445                 return True
 446             else:
 447                 return False
 448         raise NotImplementedError
 449
 450     def __ge__(self, other):
 451         return self.issuperset(other)
 452
 453     def __gt__(self, other):
 454         raise NotImplementedError
 455
 456     def union(self, *others):
 457         # return a new polyhedron with elements from the polyhedron and all
 458         # others (convex union)
 459         raise NotImplementedError
 460
 461     def __or__(self, other):
 462         return self.union(other)
 463
 464     def intersection(self, *others):
 465         # return a new polyhedron with elements common to the polyhedron and all
 466         # others
 467         # a poor man's implementation could be:
 468         # equalities = list(self.equalities)
 469         # inequalities = list(self.inequalities)
 470         # for other in others:
 471         #     equalities.extend(other.equalities)
 472         #     inequalities.extend(other.inequalities)
 473         # return self.__class__(equalities, inequalities)
 474         raise NotImplementedError
 475
 476     def __and__(self, other):
 477         return self.intersection(other)
 478
 479     def difference(self, *others):
 480         # return a new polyhedron with elements in the polyhedron that are not
 481         # in the others
 482         raise NotImplementedError
 483
 484     def __sub__(self, other):
 485         return self.difference(other)
 486
 487     def __str__(self):
 488         constraints = []
 489         for constraint in self.equalities:
 490             constraints.append('{} == 0'.format(constraint))
 491         for constraint in self.inequalities:
 492             constraints.append('{} <= 0'.format(constraint))
 493         return '{{{}}}'.format(', '.join(constraints))
 494
 495     def __repr__(self):
 496         equalities = list(self.equalities)
 497         inequalities = list(self.inequalities)
 498         return '{}(equalities={!r}, inequalities={!r})' \
 499                 ''.format(self.__class__.__name__, equalities, inequalities)
 500
 501     @classmethod
 502     def fromstring(cls, string):
 503         raise NotImplementedError
 504
 505     def to_isl(self):
 506         space = libisl.isl_space_set_alloc(_CONTEXT, 0, len(self.symbol_count()))
 507         bset = libisl.isl_basic_set_empty(libisl.isl_space_copy(space))
 508         ls = libisl.isl_local_space_from_space(libisl.isl_space_copy(space))
 509         ceq = libisl.isl_equality_alloc(libisl.isl_local_space_copy(ls))
 510         cin = libisl.isl_inequality_alloc(libisl.isl_local_space_copy(ls))
 511         d = Expression().__dict__  #write expression values to dictionary in form {'_constant': value, '_coefficients': value}
 512         '''
 513         if there are equalities/inequalities, take each constant and coefficient and add as a constraint to the basic set
 514         need to change the symbols method to a lookup table for the integer value for each letter that could be a symbol
 515         '''
 516         if self._equalities:
 517             if '_constant' in d:
 518                 value = d.get('_constant')
 519                 ceq = libisl.isl_constraint_set_constant_si(ceq, value)
 520             if '_coefficients' in d:
 521                 value_co = d.get('_coefficients')
 522                 if value_co:   #if dictionary not empty add coefficient as to constraint
 523                     ceq = libisl.isl_constraint_set_coefficient_si(ceq, libisl.isl_set_dim, self.symbols(), value_co)
 524             bset = libisl.isl_set_add_constraint(bset, ceq)
 525
 526         if self._inequalities:
 527             if '_constant' in d:
 528                 value = d.get('_constant')
 529                 cin = libisl.isl_constraint_set_constant_si(cin, value)
 530             if '_coefficients' in d:
 531                 value_co = d.get('_coefficients')
 532                 if value_co: #if dictionary not empty add coefficient as to constraint
 533                     cin = libisl.isl_constraint_set_coefficient_si(cin, libisl.isl_set_dim, self.symbols(), value_co)
 534             bset = libisl.isl_set_add_constraint(bset, cin)
 535         ip = libisl.isl_printer_to_str(_CONTEXT) #create string printer
 536         ip = libisl.isl_printer_print_set(ip, bset) #print set to printer
 537         string = libisl.isl_printer_get_str(ip)   #get string from printer
 538         print(string)
 539         return bset
 540
 541 empty = eq(1, 1)
 542
 543
 544 universe = Polyhedron()