sympycore/calculus/algebra.py

   1 #
   2 # Created January 2008 by Pearu Peterson
   3 #
   4 """ Provides Calculus class.
   5 """
   6 __docformat__ = "restructuredtext"
   7 __all__ = ['Calculus', 'I']
   8
   9 from ..core import classes, defined_functions
  10 from ..utils import TERMS, str_PRODUCT, SYMBOL, NUMBER
  11 from ..heads import APPLY, CALLABLE, TERM_COEFF, TERM_COEFF_DICT
  12 from ..heads import BASE_EXP_DICT as FACTORS
  13 from ..basealgebra import Algebra, Verbatim
  14 from ..ring import CommutativeRing
  15
  16 from ..arithmetic.numbers import normalized_fraction, mpq, mpf, mpc, mpqc, try_power
  17
  18 from ..arithmetic import mpmath, setdps
  19 from ..arithmetic.evalf import evalf
  20
  21 algebra_numbers = (int, long, mpq, mpqc, mpf, mpc)
  22 convertible_numbers = algebra_numbers + (float, complex)
  23
  24 float_one = mpf('1.0')
  25
  26 class Calculus(CommutativeRing):
  27     """ Represents an element of a symbolic algebra.
  28
  29     The set of a symbolic algebra is a set of expressions.
  30     """
  31
  32     _hash = None
  33
  34     coefftypes = algebra_numbers
  35     coefftypes_set = frozenset(coefftypes)
  36     exptypes = algebra_numbers
  37     exptypes_set = frozenset(exptypes)
  38
  39     @classmethod
  40     def get_function_algebra(cls):
  41         return classes.CalculusFunctionRing
  42
  43     def __abs__(self):
  44         head, data = self.pair
  45         if head is NUMBER:
  46             return type(self)(NUMBER, abs(data))
  47         raise NotImplementedError('abs(%r)' % (self))
  48
  49     def __call__(self, *args, **kws):
  50         head, data = self.pair
  51         if head is CALLABLE:
  52             return data(*args, **kws)
  53         raise TypeError(`self` + ' is not callable')
  54
  55     def __lt__(self, other):
  56         if type(other) is type(self):
  57             return self.data < other.data
  58         return self.data < other
  59     def __le__(self, other):
  60         if type(other) is type(self):
  61             return self.data <= other.data
  62         return self.data <= other
  63     def __gt__(self, other):
  64         if type(other) is type(self):
  65             return self.data > other.data
  66         return self.data > other
  67     def __ge__(self, other):
  68         if type(other) is type(self):
  69             return self.data >= other.data
  70         return self.data >= other
  71     def __ne__(self, other):
  72         return not (self == other)
  73
  74     def _sympy_(self):
  75         import sympy as m
  76         head, data = self.pair
  77         if head is SYMBOL:
  78             return m.Symbol(data)
  79         if head is NUMBER:
  80             return m.Number(data)
  81         if head is TERMS:
  82             return m.Add(*self.as_Add_args())
  83         if head is FACTORS:
  84             return m.Mul(*self.as_Mul_args())
  85         raise NotImplemented(`self`)
  86
  87     def as_algebra(self, cls, typeerror=True):
  88         """ Convert algebra to another algebra.
  89         """
  90         if cls is Verbatim:
  91             return self.as_verbatim()
  92         if cls is classes.Unit:
  93             return cls(NUMBER, self)
  94         if issubclass(cls, PolynomialRing):
  95             return self.as_polynom(cls)
  96         return self.as_verbatim().as_algebra(cls)
  97
  98     @classmethod
  99     def get_predefined_symbols(cls, name):
 100         if name=='I':
 101             return I
 102         if name=='D':
 103             return classes.DFactory()
 104         return getattr(defined_functions, name, None)
 105
 106     @classmethod
 107     def get_defined_function(cls, name):
 108         func = getattr(defined_functions, name, None)
 109         if func is None:
 110             func = getattr(defined_functions, name.title(), None)
 111         return func
 112
 113     @classmethod
 114     def convert_coefficient(cls, obj, typeerror=True):
 115         """ Convert obj to coefficient algebra.
 116         """
 117         if isinstance(obj, float):
 118             return mpf(obj)
 119         if isinstance(obj, complex):
 120             return mpc(obj.real, obj.imag)
 121         if isinstance(obj, algebra_numbers):
 122             return obj
 123         if isinstance(obj, cls):
 124             head, data = obj.pair
 125             if head is NUMBER:
 126                 return data
 127         if typeerror:
 128             raise TypeError('%s.convert_coefficient: failed to convert %s instance'\
 129                             ' to coefficient algebra, expected int|long object'\
 130                             % (cls.__name__, obj.__class__.__name__))
 131         else:
 132             return NotImplemented
 133
 134     @classmethod
 135     def convert_exponent(cls, obj, typeerror=True):
 136         """ Convert obj to exponent algebra.
 137         """
 138         if isinstance(obj, cls):
 139             return obj
 140         if isinstance(obj, algebra_numbers):
 141             return obj
 142         if isinstance(obj, float):
 143             return mpf(obj)
 144         if isinstance(obj, complex):
 145             return mpc(obj.real, obj.imag)
 146
 147         # parse algebra expression from string:
 148         if isinstance(obj, (str, unicode, Verbatim)):
 149             return Verbatim.convert(obj).as_algebra(cls)
 150
 151         # convert from another algebra:
 152         if isinstance(obj, Algebra):
 153             return obj.as_algebra(cls)
 154
 155         if typeerror:
 156             raise TypeError('%s.convert_exponent: failed to convert %s instance'\
 157                             ' to exponent algebra, expected int|long object'\
 158                             % (cls.__name__, obj.__class__.__name__))
 159         else:
 160             return NotImplemented
 161
 162     @classmethod
 163     def Number(cls, num, denom=None):
 164         if denom is None:
 165             return cls(NUMBER, num)
 166         return cls(NUMBER, normalized_fraction(num, denom))
 167
 168     @classmethod
 169     def Log(cls, arg, base=None):
 170         log = defined_functions.Log
 171         if base is None:
 172             return log(arg)
 173         return log(arg)/log(base)
 174
 175     @classmethod
 176     def Exp(cls, arg):
 177         return defined_functions.Exp(arg)
 178
 179     @classmethod
 180     def Mod(cls, x, y):
 181         return defined_functions.Mod(x, y)
 182
 183     def evalf(self, n=None):
 184         if n is not None:
 185             setdps(n)
 186         head, data = self.pair
 187         def evalf(cls, head, data, target):
 188             if head is NUMBER:
 189                 data = data * float_one
 190             if head is SYMBOL:
 191                 if hasattr(data, 'evalf'):
 192                     return cls(data.evalf(n))
 193             if head is APPLY:
 194                 func, args = data
 195                 name = str(func).lower()
 196                 if hasattr(mpmath, name):
 197                     f = getattr(mpmath, name)
 198                     l = []
 199                     for a in args:
 200                         h, d = a.pair
 201                         if h is NUMBER:
 202                             l.append(d)
 203                         else:
 204                             l = None
 205                             break
 206                     if l is not None:
 207                         return cls(f(*l))
 208             if head is TERM_COEFF:
 209                 term, coeff = data
 210                 coeff = coeff * float_one
 211                 return cls(head, (term, coeff))
 212             if head is TERM_COEFF_DICT:
 213                 new_data = {}
 214                 for term, coeff in data.iteritems():
 215                     new_data[term] = coeff * float_one
 216                 return cls(head, new_data)
 217             return target
 218         return head.walk(evalf, type(self), data, self)
 219
 220         if n:
 221             setdps(n)
 222         head, data = self.pair
 223         if head is NUMBER:
 224             return self.Number(data * float_one)
 225         if head is SYMBOL:
 226             r = getattr(data, 'evalf', lambda p: NotImplemented)(n)
 227             if r is not NotImplemented:
 228                 return self.Number(r)
 229             return self
 230         if head is APPLY:
 231             func = data[0]
 232             h, func_data = func.pair
 233             assert h is CALLABLE, `self, func` # todo: support symbolic functions
 234             args = data[1]
 235             assert len (args)==1,`self, args` # todo: support multivariate functions
 236             v = args[0].evalf(n)
 237             h, d = v.pair
 238             if h is NUMBER:
 239                 return self.Number(getattr(mpmath, func_data.__name__.lower())(d))
 240             else:
 241                 return type(self)(APPLY, (func, (v,)))
 242         convert = self.convert
 243         return self.func(*[convert(a).evalf(n) for a in self.args])
 244
 245     def get_direction(self):
 246         head, data = self.pair
 247         if head is NUMBER:
 248             if isinstance(data, (int, long)):
 249                 return data
 250             return getattr(data, 'get_direction', lambda : NotImplemented)()
 251         if head is TERMS:
 252             if len(data)==1:
 253                 t, c = data.items()[0]
 254                 r = t.get_direction()
 255                 if r is not NotImplemented:
 256                     return r * c
 257         if head is FACTORS:
 258             direction = 1
 259             cls = type(self)
 260             for t,c in self.data.iteritems():
 261                 d = t.get_direction()
 262                 if d is NotImplemented:
 263                     return d
 264                 if not isinstance(c, (int, long)):
 265                     return NotImplemented
 266                 d = self.Pow(cls.convert(d), c).get_direction()
 267                 if d is NotImplemented:
 268                     return d
 269                 direction *= d
 270             return direction
 271         return getattr(data, 'get_direction', lambda : NotImplemented)()
 272
 273     @property
 274     def is_bounded(self):
 275         head, data = self.pair
 276         if head is NUMBER:
 277             if isinstance(data, (int, long)):
 278                 return True
 279             return getattr(data, 'is_bounded', None)
 280         if head is SYMBOL:
 281             return getattr(data, 'is_bounded', None)
 282         if head is TERMS:
 283             for t, c in data.iteritems():
 284                 b = t.is_bounded
 285                 if not b:
 286                     return b
 287                 if isinstance(c, (int, long)):
 288                     continue
 289                 b = getattr(c, 'is_bounded', None)
 290                 if not b:
 291                     return b
 292             return True
 293         return
 294
 295     def as_polynom(self, ring_cls=None):
 296         """ Convert expression to an element of polynomial ring.
 297
 298         If the polynomial ring is not given then it will be created.
 299
 300         For example,
 301
 302           >>> x,y,z = map(Symbol,'xyz')
 303           >>> (x+2*y+3*z).as_polynom() + y*z
 304           PolynomialRing[(x, y, z), Calculus]('2*x + y*z + y + 3*z')
 305           >>> P = PolynomialRing[x,y]
 306           >>> (x+2*y+3*z).as_polynom(P) + y*z
 307           PolynomialRing[(x, y), Calculus]('x + (2 + z)*y + 3*z')
 308
 309         """
 310         if ring_cls is None:
 311             data, variables = self.to_polynomial_data()
 312             return PolynomialRing[tuple(variables)].convert(data)
 313         else:
 314             data, variables = self.to_polynomial_data(ring_cls.variables, True)
 315             return ring_cls.convert(data)
 316
 317     def __divmod__(self, other):
 318         if isinstance(other, Calculus):
 319             lhs = self.as_polynom()
 320             rhs = other.as_polynom(type(lhs))
 321             return divmod(lhs, rhs)
 322         return NotImplemented
 323
 324     def __float__(self):
 325         head, data = self.pair
 326         if head is NUMBER:
 327             return float(data)
 328         head, data = self.evalf().pair
 329         if head is NUMBER:
 330             return float(data)
 331         raise TypeError('Cannot convert %r to float' % (self))
 332
 333     def __or__(self, other):
 334         if type(other) is tuple:
 335             return self.subs(*other)
 336         return self.subs(other)
 337
 338 classes.Calculus = Calculus
 339
 340 # as a result, x<y will return Logic(LT, (x, y)):
 341 Calculus.enable_symbolic_comparison('inequality')
 342
 343 one = Calculus.Number(1)
 344 zero = Calculus.Number(0)
 345 Calculus.one = one
 346 Calculus.zero = zero
 347
 348 I = Calculus.Number(mpqc(0,1))
 349
 350 from ..polynomials.algebra import PolynomialRing, AdditiveTuple
 351 from .infinity import CalculusInfinity