Coverage for /usr/lib/python3/dist-packages/sympy/core/mul.py: 41%
1253 statements
« prev ^ index » next coverage.py v7.9.1, created at 2025-06-14 15:55 +0200
1from typing import Tuple as tTuple
2from collections import defaultdict
3from functools import cmp_to_key, reduce
4from itertools import product
5import operator
7from .sympify import sympify
8from .basic import Basic
9from .singleton import S
10from .operations import AssocOp, AssocOpDispatcher
11from .cache import cacheit
12from .logic import fuzzy_not, _fuzzy_group
13from .expr import Expr
14from .parameters import global_parameters
15from .kind import KindDispatcher
16from .traversal import bottom_up
18from sympy.utilities.iterables import sift
20# internal marker to indicate:
21# "there are still non-commutative objects -- don't forget to process them"
22class NC_Marker:
23 is_Order = False
24 is_Mul = False
25 is_Number = False
26 is_Poly = False
28 is_commutative = False
31# Key for sorting commutative args in canonical order
32_args_sortkey = cmp_to_key(Basic.compare)
33def _mulsort(args):
34 # in-place sorting of args
35 args.sort(key=_args_sortkey)
38def _unevaluated_Mul(*args):
39 """Return a well-formed unevaluated Mul: Numbers are collected and
40 put in slot 0, any arguments that are Muls will be flattened, and args
41 are sorted. Use this when args have changed but you still want to return
42 an unevaluated Mul.
44 Examples
45 ========
47 >>> from sympy.core.mul import _unevaluated_Mul as uMul
48 >>> from sympy import S, sqrt, Mul
49 >>> from sympy.abc import x
50 >>> a = uMul(*[S(3.0), x, S(2)])
51 >>> a.args[0]
52 6.00000000000000
53 >>> a.args[1]
54 x
56 Two unevaluated Muls with the same arguments will
57 always compare as equal during testing:
59 >>> m = uMul(sqrt(2), sqrt(3))
60 >>> m == uMul(sqrt(3), sqrt(2))
61 True
62 >>> u = Mul(sqrt(3), sqrt(2), evaluate=False)
63 >>> m == uMul(u)
64 True
65 >>> m == Mul(*m.args)
66 False
68 """
69 args = list(args)
70 newargs = []
71 ncargs = []
72 co = S.One
73 while args:
74 a = args.pop()
75 if a.is_Mul:
76 c, nc = a.args_cnc()
77 args.extend(c)
78 if nc:
79 ncargs.append(Mul._from_args(nc))
80 elif a.is_Number:
81 co *= a
82 else:
83 newargs.append(a)
84 _mulsort(newargs)
85 if co is not S.One:
86 newargs.insert(0, co)
87 if ncargs:
88 newargs.append(Mul._from_args(ncargs))
89 return Mul._from_args(newargs)
92class Mul(Expr, AssocOp):
93 """
94 Expression representing multiplication operation for algebraic field.
96 .. deprecated:: 1.7
98 Using arguments that aren't subclasses of :class:`~.Expr` in core
99 operators (:class:`~.Mul`, :class:`~.Add`, and :class:`~.Pow`) is
100 deprecated. See :ref:`non-expr-args-deprecated` for details.
102 Every argument of ``Mul()`` must be ``Expr``. Infix operator ``*``
103 on most scalar objects in SymPy calls this class.
105 Another use of ``Mul()`` is to represent the structure of abstract
106 multiplication so that its arguments can be substituted to return
107 different class. Refer to examples section for this.
109 ``Mul()`` evaluates the argument unless ``evaluate=False`` is passed.
110 The evaluation logic includes:
112 1. Flattening
113 ``Mul(x, Mul(y, z))`` -> ``Mul(x, y, z)``
115 2. Identity removing
116 ``Mul(x, 1, y)`` -> ``Mul(x, y)``
118 3. Exponent collecting by ``.as_base_exp()``
119 ``Mul(x, x**2)`` -> ``Pow(x, 3)``
121 4. Term sorting
122 ``Mul(y, x, 2)`` -> ``Mul(2, x, y)``
124 Since multiplication can be vector space operation, arguments may
125 have the different :obj:`sympy.core.kind.Kind()`. Kind of the
126 resulting object is automatically inferred.
128 Examples
129 ========
131 >>> from sympy import Mul
132 >>> from sympy.abc import x, y
133 >>> Mul(x, 1)
134 x
135 >>> Mul(x, x)
136 x**2
138 If ``evaluate=False`` is passed, result is not evaluated.
140 >>> Mul(1, 2, evaluate=False)
141 1*2
142 >>> Mul(x, x, evaluate=False)
143 x*x
145 ``Mul()`` also represents the general structure of multiplication
146 operation.
148 >>> from sympy import MatrixSymbol
149 >>> A = MatrixSymbol('A', 2,2)
150 >>> expr = Mul(x,y).subs({y:A})
151 >>> expr
152 x*A
153 >>> type(expr)
154 <class 'sympy.matrices.expressions.matmul.MatMul'>
156 See Also
157 ========
159 MatMul
161 """
162 __slots__ = ()
164 args: tTuple[Expr]
166 is_Mul = True
168 _args_type = Expr
169 _kind_dispatcher = KindDispatcher("Mul_kind_dispatcher", commutative=True)
171 @property
172 def kind(self):
173 arg_kinds = (a.kind for a in self.args)
174 return self._kind_dispatcher(*arg_kinds)
176 def could_extract_minus_sign(self):
177 if self == (-self):
178 return False # e.g. zoo*x == -zoo*x
179 c = self.args[0]
180 return c.is_Number and c.is_extended_negative
182 def __neg__(self):
183 c, args = self.as_coeff_mul()
184 if args[0] is not S.ComplexInfinity:
185 c = -c
186 if c is not S.One:
187 if args[0].is_Number:
188 args = list(args)
189 if c is S.NegativeOne:
190 args[0] = -args[0]
191 else:
192 args[0] *= c
193 else:
194 args = (c,) + args
195 return self._from_args(args, self.is_commutative)
197 @classmethod
198 def flatten(cls, seq):
199 """Return commutative, noncommutative and order arguments by
200 combining related terms.
202 Notes
203 =====
204 * In an expression like ``a*b*c``, Python process this through SymPy
205 as ``Mul(Mul(a, b), c)``. This can have undesirable consequences.
207 - Sometimes terms are not combined as one would like:
208 {c.f. https://github.com/sympy/sympy/issues/4596}
210 >>> from sympy import Mul, sqrt
211 >>> from sympy.abc import x, y, z
212 >>> 2*(x + 1) # this is the 2-arg Mul behavior
213 2*x + 2
214 >>> y*(x + 1)*2
215 2*y*(x + 1)
216 >>> 2*(x + 1)*y # 2-arg result will be obtained first
217 y*(2*x + 2)
218 >>> Mul(2, x + 1, y) # all 3 args simultaneously processed
219 2*y*(x + 1)
220 >>> 2*((x + 1)*y) # parentheses can control this behavior
221 2*y*(x + 1)
223 Powers with compound bases may not find a single base to
224 combine with unless all arguments are processed at once.
225 Post-processing may be necessary in such cases.
226 {c.f. https://github.com/sympy/sympy/issues/5728}
228 >>> a = sqrt(x*sqrt(y))
229 >>> a**3
230 (x*sqrt(y))**(3/2)
231 >>> Mul(a,a,a)
232 (x*sqrt(y))**(3/2)
233 >>> a*a*a
234 x*sqrt(y)*sqrt(x*sqrt(y))
235 >>> _.subs(a.base, z).subs(z, a.base)
236 (x*sqrt(y))**(3/2)
238 - If more than two terms are being multiplied then all the
239 previous terms will be re-processed for each new argument.
240 So if each of ``a``, ``b`` and ``c`` were :class:`Mul`
241 expression, then ``a*b*c`` (or building up the product
242 with ``*=``) will process all the arguments of ``a`` and
243 ``b`` twice: once when ``a*b`` is computed and again when
244 ``c`` is multiplied.
246 Using ``Mul(a, b, c)`` will process all arguments once.
248 * The results of Mul are cached according to arguments, so flatten
249 will only be called once for ``Mul(a, b, c)``. If you can
250 structure a calculation so the arguments are most likely to be
251 repeats then this can save time in computing the answer. For
252 example, say you had a Mul, M, that you wished to divide by ``d[i]``
253 and multiply by ``n[i]`` and you suspect there are many repeats
254 in ``n``. It would be better to compute ``M*n[i]/d[i]`` rather
255 than ``M/d[i]*n[i]`` since every time n[i] is a repeat, the
256 product, ``M*n[i]`` will be returned without flattening -- the
257 cached value will be returned. If you divide by the ``d[i]``
258 first (and those are more unique than the ``n[i]``) then that will
259 create a new Mul, ``M/d[i]`` the args of which will be traversed
260 again when it is multiplied by ``n[i]``.
262 {c.f. https://github.com/sympy/sympy/issues/5706}
264 This consideration is moot if the cache is turned off.
266 NB
267 --
268 The validity of the above notes depends on the implementation
269 details of Mul and flatten which may change at any time. Therefore,
270 you should only consider them when your code is highly performance
271 sensitive.
273 Removal of 1 from the sequence is already handled by AssocOp.__new__.
274 """
276 from sympy.calculus.accumulationbounds import AccumBounds
277 from sympy.matrices.expressions import MatrixExpr
278 rv = None
279 if len(seq) == 2:
280 a, b = seq
281 if b.is_Rational:
282 a, b = b, a
283 seq = [a, b]
284 assert a is not S.One
285 if not a.is_zero and a.is_Rational:
286 r, b = b.as_coeff_Mul()
287 if b.is_Add:
288 if r is not S.One: # 2-arg hack
289 # leave the Mul as a Mul?
290 ar = a*r
291 if ar is S.One:
292 arb = b
293 else:
294 arb = cls(a*r, b, evaluate=False)
295 rv = [arb], [], None
296 elif global_parameters.distribute and b.is_commutative:
297 newb = Add(*[_keep_coeff(a, bi) for bi in b.args])
298 rv = [newb], [], None
299 if rv:
300 return rv
302 # apply associativity, separate commutative part of seq
303 c_part = [] # out: commutative factors
304 nc_part = [] # out: non-commutative factors
306 nc_seq = []
308 coeff = S.One # standalone term
309 # e.g. 3 * ...
311 c_powers = [] # (base,exp) n
312 # e.g. (x,n) for x
314 num_exp = [] # (num-base, exp) y
315 # e.g. (3, y) for ... * 3 * ...
317 neg1e = S.Zero # exponent on -1 extracted from Number-based Pow and I
319 pnum_rat = {} # (num-base, Rat-exp) 1/2
320 # e.g. (3, 1/2) for ... * 3 * ...
322 order_symbols = None
324 # --- PART 1 ---
325 #
326 # "collect powers and coeff":
327 #
328 # o coeff
329 # o c_powers
330 # o num_exp
331 # o neg1e
332 # o pnum_rat
333 #
334 # NOTE: this is optimized for all-objects-are-commutative case
335 for o in seq:
336 # O(x)
337 if o.is_Order:
338 o, order_symbols = o.as_expr_variables(order_symbols)
340 # Mul([...])
341 if o.is_Mul:
342 if o.is_commutative:
343 seq.extend(o.args) # XXX zerocopy?
345 else:
346 # NCMul can have commutative parts as well
347 for q in o.args:
348 if q.is_commutative:
349 seq.append(q)
350 else:
351 nc_seq.append(q)
353 # append non-commutative marker, so we don't forget to
354 # process scheduled non-commutative objects
355 seq.append(NC_Marker)
357 continue
359 # 3
360 elif o.is_Number:
361 if o is S.NaN or coeff is S.ComplexInfinity and o.is_zero:
362 # we know for sure the result will be nan
363 return [S.NaN], [], None
364 elif coeff.is_Number or isinstance(coeff, AccumBounds): # it could be zoo
365 coeff *= o
366 if coeff is S.NaN:
367 # we know for sure the result will be nan
368 return [S.NaN], [], None
369 continue
371 elif isinstance(o, AccumBounds):
372 coeff = o.__mul__(coeff)
373 continue
375 elif o is S.ComplexInfinity:
376 if not coeff:
377 # 0 * zoo = NaN
378 return [S.NaN], [], None
379 coeff = S.ComplexInfinity
380 continue
382 elif o is S.ImaginaryUnit:
383 neg1e += S.Half
384 continue
386 elif o.is_commutative:
387 # e
388 # o = b
389 b, e = o.as_base_exp()
391 # y
392 # 3
393 if o.is_Pow:
394 if b.is_Number:
396 # get all the factors with numeric base so they can be
397 # combined below, but don't combine negatives unless
398 # the exponent is an integer
399 if e.is_Rational:
400 if e.is_Integer:
401 coeff *= Pow(b, e) # it is an unevaluated power
402 continue
403 elif e.is_negative: # also a sign of an unevaluated power
404 seq.append(Pow(b, e))
405 continue
406 elif b.is_negative:
407 neg1e += e
408 b = -b
409 if b is not S.One:
410 pnum_rat.setdefault(b, []).append(e)
411 continue
412 elif b.is_positive or e.is_integer:
413 num_exp.append((b, e))
414 continue
416 c_powers.append((b, e))
418 # NON-COMMUTATIVE
419 # TODO: Make non-commutative exponents not combine automatically
420 else:
421 if o is not NC_Marker:
422 nc_seq.append(o)
424 # process nc_seq (if any)
425 while nc_seq:
426 o = nc_seq.pop(0)
427 if not nc_part:
428 nc_part.append(o)
429 continue
431 # b c b+c
432 # try to combine last terms: a * a -> a
433 o1 = nc_part.pop()
434 b1, e1 = o1.as_base_exp()
435 b2, e2 = o.as_base_exp()
436 new_exp = e1 + e2
437 # Only allow powers to combine if the new exponent is
438 # not an Add. This allow things like a**2*b**3 == a**5
439 # if a.is_commutative == False, but prohibits
440 # a**x*a**y and x**a*x**b from combining (x,y commute).
441 if b1 == b2 and (not new_exp.is_Add):
442 o12 = b1 ** new_exp
444 # now o12 could be a commutative object
445 if o12.is_commutative:
446 seq.append(o12)
447 continue
448 else:
449 nc_seq.insert(0, o12)
451 else:
452 nc_part.extend([o1, o])
454 # We do want a combined exponent if it would not be an Add, such as
455 # y 2y 3y
456 # x * x -> x
457 # We determine if two exponents have the same term by using
458 # as_coeff_Mul.
459 #
460 # Unfortunately, this isn't smart enough to consider combining into
461 # exponents that might already be adds, so things like:
462 # z - y y
463 # x * x will be left alone. This is because checking every possible
464 # combination can slow things down.
466 # gather exponents of common bases...
467 def _gather(c_powers):
468 common_b = {} # b:e
469 for b, e in c_powers:
470 co = e.as_coeff_Mul()
471 common_b.setdefault(b, {}).setdefault(
472 co[1], []).append(co[0])
473 for b, d in common_b.items():
474 for di, li in d.items():
475 d[di] = Add(*li)
476 new_c_powers = []
477 for b, e in common_b.items():
478 new_c_powers.extend([(b, c*t) for t, c in e.items()])
479 return new_c_powers
481 # in c_powers
482 c_powers = _gather(c_powers)
484 # and in num_exp
485 num_exp = _gather(num_exp)
487 # --- PART 2 ---
488 #
489 # o process collected powers (x**0 -> 1; x**1 -> x; otherwise Pow)
490 # o combine collected powers (2**x * 3**x -> 6**x)
491 # with numeric base
493 # ................................
494 # now we have:
495 # - coeff:
496 # - c_powers: (b, e)
497 # - num_exp: (2, e)
498 # - pnum_rat: {(1/3, [1/3, 2/3, 1/4])}
500 # 0 1
501 # x -> 1 x -> x
503 # this should only need to run twice; if it fails because
504 # it needs to be run more times, perhaps this should be
505 # changed to a "while True" loop -- the only reason it
506 # isn't such now is to allow a less-than-perfect result to
507 # be obtained rather than raising an error or entering an
508 # infinite loop
509 for i in range(2):
510 new_c_powers = []
511 changed = False
512 for b, e in c_powers:
513 if e.is_zero:
514 # canceling out infinities yields NaN
515 if (b.is_Add or b.is_Mul) and any(infty in b.args
516 for infty in (S.ComplexInfinity, S.Infinity,
517 S.NegativeInfinity)):
518 return [S.NaN], [], None
519 continue
520 if e is S.One:
521 if b.is_Number:
522 coeff *= b
523 continue
524 p = b
525 if e is not S.One:
526 p = Pow(b, e)
527 # check to make sure that the base doesn't change
528 # after exponentiation; to allow for unevaluated
529 # Pow, we only do so if b is not already a Pow
530 if p.is_Pow and not b.is_Pow:
531 bi = b
532 b, e = p.as_base_exp()
533 if b != bi:
534 changed = True
535 c_part.append(p)
536 new_c_powers.append((b, e))
537 # there might have been a change, but unless the base
538 # matches some other base, there is nothing to do
539 if changed and len({
540 b for b, e in new_c_powers}) != len(new_c_powers):
541 # start over again
542 c_part = []
543 c_powers = _gather(new_c_powers)
544 else:
545 break
547 # x x x
548 # 2 * 3 -> 6
549 inv_exp_dict = {} # exp:Mul(num-bases) x x
550 # e.g. x:6 for ... * 2 * 3 * ...
551 for b, e in num_exp:
552 inv_exp_dict.setdefault(e, []).append(b)
553 for e, b in inv_exp_dict.items():
554 inv_exp_dict[e] = cls(*b)
555 c_part.extend([Pow(b, e) for e, b in inv_exp_dict.items() if e])
557 # b, e -> e' = sum(e), b
558 # {(1/5, [1/3]), (1/2, [1/12, 1/4]} -> {(1/3, [1/5, 1/2])}
559 comb_e = {}
560 for b, e in pnum_rat.items():
561 comb_e.setdefault(Add(*e), []).append(b)
562 del pnum_rat
563 # process them, reducing exponents to values less than 1
564 # and updating coeff if necessary else adding them to
565 # num_rat for further processing
566 num_rat = []
567 for e, b in comb_e.items():
568 b = cls(*b)
569 if e.q == 1:
570 coeff *= Pow(b, e)
571 continue
572 if e.p > e.q:
573 e_i, ep = divmod(e.p, e.q)
574 coeff *= Pow(b, e_i)
575 e = Rational(ep, e.q)
576 num_rat.append((b, e))
577 del comb_e
579 # extract gcd of bases in num_rat
580 # 2**(1/3)*6**(1/4) -> 2**(1/3+1/4)*3**(1/4)
581 pnew = defaultdict(list)
582 i = 0 # steps through num_rat which may grow
583 while i < len(num_rat):
584 bi, ei = num_rat[i]
585 grow = []
586 for j in range(i + 1, len(num_rat)):
587 bj, ej = num_rat[j]
588 g = bi.gcd(bj)
589 if g is not S.One:
590 # 4**r1*6**r2 -> 2**(r1+r2) * 2**r1 * 3**r2
591 # this might have a gcd with something else
592 e = ei + ej
593 if e.q == 1:
594 coeff *= Pow(g, e)
595 else:
596 if e.p > e.q:
597 e_i, ep = divmod(e.p, e.q) # change e in place
598 coeff *= Pow(g, e_i)
599 e = Rational(ep, e.q)
600 grow.append((g, e))
601 # update the jth item
602 num_rat[j] = (bj/g, ej)
603 # update bi that we are checking with
604 bi = bi/g
605 if bi is S.One:
606 break
607 if bi is not S.One:
608 obj = Pow(bi, ei)
609 if obj.is_Number:
610 coeff *= obj
611 else:
612 # changes like sqrt(12) -> 2*sqrt(3)
613 for obj in Mul.make_args(obj):
614 if obj.is_Number:
615 coeff *= obj
616 else:
617 assert obj.is_Pow
618 bi, ei = obj.args
619 pnew[ei].append(bi)
621 num_rat.extend(grow)
622 i += 1
624 # combine bases of the new powers
625 for e, b in pnew.items():
626 pnew[e] = cls(*b)
628 # handle -1 and I
629 if neg1e:
630 # treat I as (-1)**(1/2) and compute -1's total exponent
631 p, q = neg1e.as_numer_denom()
632 # if the integer part is odd, extract -1
633 n, p = divmod(p, q)
634 if n % 2:
635 coeff = -coeff
636 # if it's a multiple of 1/2 extract I
637 if q == 2:
638 c_part.append(S.ImaginaryUnit)
639 elif p:
640 # see if there is any positive base this power of
641 # -1 can join
642 neg1e = Rational(p, q)
643 for e, b in pnew.items():
644 if e == neg1e and b.is_positive:
645 pnew[e] = -b
646 break
647 else:
648 # keep it separate; we've already evaluated it as
649 # much as possible so evaluate=False
650 c_part.append(Pow(S.NegativeOne, neg1e, evaluate=False))
652 # add all the pnew powers
653 c_part.extend([Pow(b, e) for e, b in pnew.items()])
655 # oo, -oo
656 if coeff in (S.Infinity, S.NegativeInfinity):
657 def _handle_for_oo(c_part, coeff_sign):
658 new_c_part = []
659 for t in c_part:
660 if t.is_extended_positive:
661 continue
662 if t.is_extended_negative:
663 coeff_sign *= -1
664 continue
665 new_c_part.append(t)
666 return new_c_part, coeff_sign
667 c_part, coeff_sign = _handle_for_oo(c_part, 1)
668 nc_part, coeff_sign = _handle_for_oo(nc_part, coeff_sign)
669 coeff *= coeff_sign
671 # zoo
672 if coeff is S.ComplexInfinity:
673 # zoo might be
674 # infinite_real + bounded_im
675 # bounded_real + infinite_im
676 # infinite_real + infinite_im
677 # and non-zero real or imaginary will not change that status.
678 c_part = [c for c in c_part if not (fuzzy_not(c.is_zero) and
679 c.is_extended_real is not None)]
680 nc_part = [c for c in nc_part if not (fuzzy_not(c.is_zero) and
681 c.is_extended_real is not None)]
683 # 0
684 elif coeff.is_zero:
685 # we know for sure the result will be 0 except the multiplicand
686 # is infinity or a matrix
687 if any(isinstance(c, MatrixExpr) for c in nc_part):
688 return [coeff], nc_part, order_symbols
689 if any(c.is_finite == False for c in c_part):
690 return [S.NaN], [], order_symbols
691 return [coeff], [], order_symbols
693 # check for straggling Numbers that were produced
694 _new = []
695 for i in c_part:
696 if i.is_Number:
697 coeff *= i
698 else:
699 _new.append(i)
700 c_part = _new
702 # order commutative part canonically
703 _mulsort(c_part)
705 # current code expects coeff to be always in slot-0
706 if coeff is not S.One:
707 c_part.insert(0, coeff)
709 # we are done
710 if (global_parameters.distribute and not nc_part and len(c_part) == 2 and
711 c_part[0].is_Number and c_part[0].is_finite and c_part[1].is_Add):
712 # 2*(1+a) -> 2 + 2 * a
713 coeff = c_part[0]
714 c_part = [Add(*[coeff*f for f in c_part[1].args])]
716 return c_part, nc_part, order_symbols
718 def _eval_power(self, e):
720 # don't break up NC terms: (A*B)**3 != A**3*B**3, it is A*B*A*B*A*B
721 cargs, nc = self.args_cnc(split_1=False)
723 if e.is_Integer:
724 return Mul(*[Pow(b, e, evaluate=False) for b in cargs]) * \
725 Pow(Mul._from_args(nc), e, evaluate=False)
726 if e.is_Rational and e.q == 2:
727 if self.is_imaginary:
728 a = self.as_real_imag()[1]
729 if a.is_Rational:
730 from .power import integer_nthroot
731 n, d = abs(a/2).as_numer_denom()
732 n, t = integer_nthroot(n, 2)
733 if t:
734 d, t = integer_nthroot(d, 2)
735 if t:
736 from sympy.functions.elementary.complexes import sign
737 r = sympify(n)/d
738 return _unevaluated_Mul(r**e.p, (1 + sign(a)*S.ImaginaryUnit)**e.p)
740 p = Pow(self, e, evaluate=False)
742 if e.is_Rational or e.is_Float:
743 return p._eval_expand_power_base()
745 return p
747 @classmethod
748 def class_key(cls):
749 return 3, 0, cls.__name__
751 def _eval_evalf(self, prec):
752 c, m = self.as_coeff_Mul()
753 if c is S.NegativeOne:
754 if m.is_Mul:
755 rv = -AssocOp._eval_evalf(m, prec)
756 else:
757 mnew = m._eval_evalf(prec)
758 if mnew is not None:
759 m = mnew
760 rv = -m
761 else:
762 rv = AssocOp._eval_evalf(self, prec)
763 if rv.is_number:
764 return rv.expand()
765 return rv
767 @property
768 def _mpc_(self):
769 """
770 Convert self to an mpmath mpc if possible
771 """
772 from .numbers import Float
773 im_part, imag_unit = self.as_coeff_Mul()
774 if imag_unit is not S.ImaginaryUnit:
775 # ValueError may seem more reasonable but since it's a @property,
776 # we need to use AttributeError to keep from confusing things like
777 # hasattr.
778 raise AttributeError("Cannot convert Mul to mpc. Must be of the form Number*I")
780 return (Float(0)._mpf_, Float(im_part)._mpf_)
782 @cacheit
783 def as_two_terms(self):
784 """Return head and tail of self.
786 This is the most efficient way to get the head and tail of an
787 expression.
789 - if you want only the head, use self.args[0];
790 - if you want to process the arguments of the tail then use
791 self.as_coef_mul() which gives the head and a tuple containing
792 the arguments of the tail when treated as a Mul.
793 - if you want the coefficient when self is treated as an Add
794 then use self.as_coeff_add()[0]
796 Examples
797 ========
799 >>> from sympy.abc import x, y
800 >>> (3*x*y).as_two_terms()
801 (3, x*y)
802 """
803 args = self.args
805 if len(args) == 1:
806 return S.One, self
807 elif len(args) == 2:
808 return args
810 else:
811 return args[0], self._new_rawargs(*args[1:])
813 @cacheit
814 def as_coeff_mul(self, *deps, rational=True, **kwargs):
815 if deps:
816 l1, l2 = sift(self.args, lambda x: x.has(*deps), binary=True)
817 return self._new_rawargs(*l2), tuple(l1)
818 args = self.args
819 if args[0].is_Number:
820 if not rational or args[0].is_Rational:
821 return args[0], args[1:]
822 elif args[0].is_extended_negative:
823 return S.NegativeOne, (-args[0],) + args[1:]
824 return S.One, args
826 def as_coeff_Mul(self, rational=False):
827 """
828 Efficiently extract the coefficient of a product.
829 """
830 coeff, args = self.args[0], self.args[1:]
832 if coeff.is_Number:
833 if not rational or coeff.is_Rational:
834 if len(args) == 1:
835 return coeff, args[0]
836 else:
837 return coeff, self._new_rawargs(*args)
838 elif coeff.is_extended_negative:
839 return S.NegativeOne, self._new_rawargs(*((-coeff,) + args))
840 return S.One, self
842 def as_real_imag(self, deep=True, **hints):
843 from sympy.functions.elementary.complexes import Abs, im, re
844 other = []
845 coeffr = []
846 coeffi = []
847 addterms = S.One
848 for a in self.args:
849 r, i = a.as_real_imag()
850 if i.is_zero:
851 coeffr.append(r)
852 elif r.is_zero:
853 coeffi.append(i*S.ImaginaryUnit)
854 elif a.is_commutative:
855 aconj = a.conjugate() if other else None
856 # search for complex conjugate pairs:
857 for i, x in enumerate(other):
858 if x == aconj:
859 coeffr.append(Abs(x)**2)
860 del other[i]
861 break
862 else:
863 if a.is_Add:
864 addterms *= a
865 else:
866 other.append(a)
867 else:
868 other.append(a)
869 m = self.func(*other)
870 if hints.get('ignore') == m:
871 return
872 if len(coeffi) % 2:
873 imco = im(coeffi.pop(0))
874 # all other pairs make a real factor; they will be
875 # put into reco below
876 else:
877 imco = S.Zero
878 reco = self.func(*(coeffr + coeffi))
879 r, i = (reco*re(m), reco*im(m))
880 if addterms == 1:
881 if m == 1:
882 if imco.is_zero:
883 return (reco, S.Zero)
884 else:
885 return (S.Zero, reco*imco)
886 if imco is S.Zero:
887 return (r, i)
888 return (-imco*i, imco*r)
889 from .function import expand_mul
890 addre, addim = expand_mul(addterms, deep=False).as_real_imag()
891 if imco is S.Zero:
892 return (r*addre - i*addim, i*addre + r*addim)
893 else:
894 r, i = -imco*i, imco*r
895 return (r*addre - i*addim, r*addim + i*addre)
897 @staticmethod
898 def _expandsums(sums):
899 """
900 Helper function for _eval_expand_mul.
902 sums must be a list of instances of Basic.
903 """
905 L = len(sums)
906 if L == 1:
907 return sums[0].args
908 terms = []
909 left = Mul._expandsums(sums[:L//2])
910 right = Mul._expandsums(sums[L//2:])
912 terms = [Mul(a, b) for a in left for b in right]
913 added = Add(*terms)
914 return Add.make_args(added) # it may have collapsed down to one term
916 def _eval_expand_mul(self, **hints):
917 from sympy.simplify.radsimp import fraction
919 # Handle things like 1/(x*(x + 1)), which are automatically converted
920 # to 1/x*1/(x + 1)
921 expr = self
922 n, d = fraction(expr)
923 if d.is_Mul:
924 n, d = [i._eval_expand_mul(**hints) if i.is_Mul else i
925 for i in (n, d)]
926 expr = n/d
927 if not expr.is_Mul:
928 return expr
930 plain, sums, rewrite = [], [], False
931 for factor in expr.args:
932 if factor.is_Add:
933 sums.append(factor)
934 rewrite = True
935 else:
936 if factor.is_commutative:
937 plain.append(factor)
938 else:
939 sums.append(Basic(factor)) # Wrapper
941 if not rewrite:
942 return expr
943 else:
944 plain = self.func(*plain)
945 if sums:
946 deep = hints.get("deep", False)
947 terms = self.func._expandsums(sums)
948 args = []
949 for term in terms:
950 t = self.func(plain, term)
951 if t.is_Mul and any(a.is_Add for a in t.args) and deep:
952 t = t._eval_expand_mul()
953 args.append(t)
954 return Add(*args)
955 else:
956 return plain
958 @cacheit
959 def _eval_derivative(self, s):
960 args = list(self.args)
961 terms = []
962 for i in range(len(args)):
963 d = args[i].diff(s)
964 if d:
965 # Note: reduce is used in step of Mul as Mul is unable to
966 # handle subtypes and operation priority:
967 terms.append(reduce(lambda x, y: x*y, (args[:i] + [d] + args[i + 1:]), S.One))
968 return Add.fromiter(terms)
970 @cacheit
971 def _eval_derivative_n_times(self, s, n):
972 from .function import AppliedUndef
973 from .symbol import Symbol, symbols, Dummy
974 if not isinstance(s, (AppliedUndef, Symbol)):
975 # other types of s may not be well behaved, e.g.
976 # (cos(x)*sin(y)).diff([[x, y, z]])
977 return super()._eval_derivative_n_times(s, n)
978 from .numbers import Integer
979 args = self.args
980 m = len(args)
981 if isinstance(n, (int, Integer)):
982 # https://en.wikipedia.org/wiki/General_Leibniz_rule#More_than_two_factors
983 terms = []
984 from sympy.ntheory.multinomial import multinomial_coefficients_iterator
985 for kvals, c in multinomial_coefficients_iterator(m, n):
986 p = Mul(*[arg.diff((s, k)) for k, arg in zip(kvals, args)])
987 terms.append(c * p)
988 return Add(*terms)
989 from sympy.concrete.summations import Sum
990 from sympy.functions.combinatorial.factorials import factorial
991 from sympy.functions.elementary.miscellaneous import Max
992 kvals = symbols("k1:%i" % m, cls=Dummy)
993 klast = n - sum(kvals)
994 nfact = factorial(n)
995 e, l = (# better to use the multinomial?
996 nfact/prod(map(factorial, kvals))/factorial(klast)*\
997 Mul(*[args[t].diff((s, kvals[t])) for t in range(m-1)])*\
998 args[-1].diff((s, Max(0, klast))),
999 [(k, 0, n) for k in kvals])
1000 return Sum(e, *l)
1002 def _eval_difference_delta(self, n, step):
1003 from sympy.series.limitseq import difference_delta as dd
1004 arg0 = self.args[0]
1005 rest = Mul(*self.args[1:])
1006 return (arg0.subs(n, n + step) * dd(rest, n, step) + dd(arg0, n, step) *
1007 rest)
1009 def _matches_simple(self, expr, repl_dict):
1010 # handle (w*3).matches('x*5') -> {w: x*5/3}
1011 coeff, terms = self.as_coeff_Mul()
1012 terms = Mul.make_args(terms)
1013 if len(terms) == 1:
1014 newexpr = self.__class__._combine_inverse(expr, coeff)
1015 return terms[0].matches(newexpr, repl_dict)
1016 return
1018 def matches(self, expr, repl_dict=None, old=False):
1019 expr = sympify(expr)
1020 if self.is_commutative and expr.is_commutative:
1021 return self._matches_commutative(expr, repl_dict, old)
1022 elif self.is_commutative is not expr.is_commutative:
1023 return None
1025 # Proceed only if both both expressions are non-commutative
1026 c1, nc1 = self.args_cnc()
1027 c2, nc2 = expr.args_cnc()
1028 c1, c2 = [c or [1] for c in [c1, c2]]
1030 # TODO: Should these be self.func?
1031 comm_mul_self = Mul(*c1)
1032 comm_mul_expr = Mul(*c2)
1034 repl_dict = comm_mul_self.matches(comm_mul_expr, repl_dict, old)
1036 # If the commutative arguments didn't match and aren't equal, then
1037 # then the expression as a whole doesn't match
1038 if not repl_dict and c1 != c2:
1039 return None
1041 # Now match the non-commutative arguments, expanding powers to
1042 # multiplications
1043 nc1 = Mul._matches_expand_pows(nc1)
1044 nc2 = Mul._matches_expand_pows(nc2)
1046 repl_dict = Mul._matches_noncomm(nc1, nc2, repl_dict)
1048 return repl_dict or None
1050 @staticmethod
1051 def _matches_expand_pows(arg_list):
1052 new_args = []
1053 for arg in arg_list:
1054 if arg.is_Pow and arg.exp > 0:
1055 new_args.extend([arg.base] * arg.exp)
1056 else:
1057 new_args.append(arg)
1058 return new_args
1060 @staticmethod
1061 def _matches_noncomm(nodes, targets, repl_dict=None):
1062 """Non-commutative multiplication matcher.
1064 `nodes` is a list of symbols within the matcher multiplication
1065 expression, while `targets` is a list of arguments in the
1066 multiplication expression being matched against.
1067 """
1068 if repl_dict is None:
1069 repl_dict = {}
1070 else:
1071 repl_dict = repl_dict.copy()
1073 # List of possible future states to be considered
1074 agenda = []
1075 # The current matching state, storing index in nodes and targets
1076 state = (0, 0)
1077 node_ind, target_ind = state
1078 # Mapping between wildcard indices and the index ranges they match
1079 wildcard_dict = {}
1081 while target_ind < len(targets) and node_ind < len(nodes):
1082 node = nodes[node_ind]
1084 if node.is_Wild:
1085 Mul._matches_add_wildcard(wildcard_dict, state)
1087 states_matches = Mul._matches_new_states(wildcard_dict, state,
1088 nodes, targets)
1089 if states_matches:
1090 new_states, new_matches = states_matches
1091 agenda.extend(new_states)
1092 if new_matches:
1093 for match in new_matches:
1094 repl_dict[match] = new_matches[match]
1095 if not agenda:
1096 return None
1097 else:
1098 state = agenda.pop()
1099 node_ind, target_ind = state
1101 return repl_dict
1103 @staticmethod
1104 def _matches_add_wildcard(dictionary, state):
1105 node_ind, target_ind = state
1106 if node_ind in dictionary:
1107 begin, end = dictionary[node_ind]
1108 dictionary[node_ind] = (begin, target_ind)
1109 else:
1110 dictionary[node_ind] = (target_ind, target_ind)
1112 @staticmethod
1113 def _matches_new_states(dictionary, state, nodes, targets):
1114 node_ind, target_ind = state
1115 node = nodes[node_ind]
1116 target = targets[target_ind]
1118 # Don't advance at all if we've exhausted the targets but not the nodes
1119 if target_ind >= len(targets) - 1 and node_ind < len(nodes) - 1:
1120 return None
1122 if node.is_Wild:
1123 match_attempt = Mul._matches_match_wilds(dictionary, node_ind,
1124 nodes, targets)
1125 if match_attempt:
1126 # If the same node has been matched before, don't return
1127 # anything if the current match is diverging from the previous
1128 # match
1129 other_node_inds = Mul._matches_get_other_nodes(dictionary,
1130 nodes, node_ind)
1131 for ind in other_node_inds:
1132 other_begin, other_end = dictionary[ind]
1133 curr_begin, curr_end = dictionary[node_ind]
1135 other_targets = targets[other_begin:other_end + 1]
1136 current_targets = targets[curr_begin:curr_end + 1]
1138 for curr, other in zip(current_targets, other_targets):
1139 if curr != other:
1140 return None
1142 # A wildcard node can match more than one target, so only the
1143 # target index is advanced
1144 new_state = [(node_ind, target_ind + 1)]
1145 # Only move on to the next node if there is one
1146 if node_ind < len(nodes) - 1:
1147 new_state.append((node_ind + 1, target_ind + 1))
1148 return new_state, match_attempt
1149 else:
1150 # If we're not at a wildcard, then make sure we haven't exhausted
1151 # nodes but not targets, since in this case one node can only match
1152 # one target
1153 if node_ind >= len(nodes) - 1 and target_ind < len(targets) - 1:
1154 return None
1156 match_attempt = node.matches(target)
1158 if match_attempt:
1159 return [(node_ind + 1, target_ind + 1)], match_attempt
1160 elif node == target:
1161 return [(node_ind + 1, target_ind + 1)], None
1162 else:
1163 return None
1165 @staticmethod
1166 def _matches_match_wilds(dictionary, wildcard_ind, nodes, targets):
1167 """Determine matches of a wildcard with sub-expression in `target`."""
1168 wildcard = nodes[wildcard_ind]
1169 begin, end = dictionary[wildcard_ind]
1170 terms = targets[begin:end + 1]
1171 # TODO: Should this be self.func?
1172 mult = Mul(*terms) if len(terms) > 1 else terms[0]
1173 return wildcard.matches(mult)
1175 @staticmethod
1176 def _matches_get_other_nodes(dictionary, nodes, node_ind):
1177 """Find other wildcards that may have already been matched."""
1178 ind_node = nodes[node_ind]
1179 return [ind for ind in dictionary if nodes[ind] == ind_node]
1181 @staticmethod
1182 def _combine_inverse(lhs, rhs):
1183 """
1184 Returns lhs/rhs, but treats arguments like symbols, so things
1185 like oo/oo return 1 (instead of a nan) and ``I`` behaves like
1186 a symbol instead of sqrt(-1).
1187 """
1188 from sympy.simplify.simplify import signsimp
1189 from .symbol import Dummy
1190 if lhs == rhs:
1191 return S.One
1193 def check(l, r):
1194 if l.is_Float and r.is_comparable:
1195 # if both objects are added to 0 they will share the same "normalization"
1196 # and are more likely to compare the same. Since Add(foo, 0) will not allow
1197 # the 0 to pass, we use __add__ directly.
1198 return l.__add__(0) == r.evalf().__add__(0)
1199 return False
1200 if check(lhs, rhs) or check(rhs, lhs):
1201 return S.One
1202 if any(i.is_Pow or i.is_Mul for i in (lhs, rhs)):
1203 # gruntz and limit wants a literal I to not combine
1204 # with a power of -1
1205 d = Dummy('I')
1206 _i = {S.ImaginaryUnit: d}
1207 i_ = {d: S.ImaginaryUnit}
1208 a = lhs.xreplace(_i).as_powers_dict()
1209 b = rhs.xreplace(_i).as_powers_dict()
1210 blen = len(b)
1211 for bi in tuple(b.keys()):
1212 if bi in a:
1213 a[bi] -= b.pop(bi)
1214 if not a[bi]:
1215 a.pop(bi)
1216 if len(b) != blen:
1217 lhs = Mul(*[k**v for k, v in a.items()]).xreplace(i_)
1218 rhs = Mul(*[k**v for k, v in b.items()]).xreplace(i_)
1219 rv = lhs/rhs
1220 srv = signsimp(rv)
1221 return srv if srv.is_Number else rv
1223 def as_powers_dict(self):
1224 d = defaultdict(int)
1225 for term in self.args:
1226 for b, e in term.as_powers_dict().items():
1227 d[b] += e
1228 return d
1230 def as_numer_denom(self):
1231 # don't use _from_args to rebuild the numerators and denominators
1232 # as the order is not guaranteed to be the same once they have
1233 # been separated from each other
1234 numers, denoms = list(zip(*[f.as_numer_denom() for f in self.args]))
1235 return self.func(*numers), self.func(*denoms)
1237 def as_base_exp(self):
1238 e1 = None
1239 bases = []
1240 nc = 0
1241 for m in self.args:
1242 b, e = m.as_base_exp()
1243 if not b.is_commutative:
1244 nc += 1
1245 if e1 is None:
1246 e1 = e
1247 elif e != e1 or nc > 1:
1248 return self, S.One
1249 bases.append(b)
1250 return self.func(*bases), e1
1252 def _eval_is_polynomial(self, syms):
1253 return all(term._eval_is_polynomial(syms) for term in self.args)
1255 def _eval_is_rational_function(self, syms):
1256 return all(term._eval_is_rational_function(syms) for term in self.args)
1258 def _eval_is_meromorphic(self, x, a):
1259 return _fuzzy_group((arg.is_meromorphic(x, a) for arg in self.args),
1260 quick_exit=True)
1262 def _eval_is_algebraic_expr(self, syms):
1263 return all(term._eval_is_algebraic_expr(syms) for term in self.args)
1265 _eval_is_commutative = lambda self: _fuzzy_group(
1266 a.is_commutative for a in self.args)
1268 def _eval_is_complex(self):
1269 comp = _fuzzy_group(a.is_complex for a in self.args)
1270 if comp is False:
1271 if any(a.is_infinite for a in self.args):
1272 if any(a.is_zero is not False for a in self.args):
1273 return None
1274 return False
1275 return comp
1277 def _eval_is_zero_infinite_helper(self):
1278 #
1279 # Helper used by _eval_is_zero and _eval_is_infinite.
1280 #
1281 # Three-valued logic is tricky so let us reason this carefully. It
1282 # would be nice to say that we just check is_zero/is_infinite in all
1283 # args but we need to be careful about the case that one arg is zero
1284 # and another is infinite like Mul(0, oo) or more importantly a case
1285 # where it is not known if the arguments are zero or infinite like
1286 # Mul(y, 1/x). If either y or x could be zero then there is a
1287 # *possibility* that we have Mul(0, oo) which should give None for both
1288 # is_zero and is_infinite.
1289 #
1290 # We keep track of whether we have seen a zero or infinity but we also
1291 # need to keep track of whether we have *possibly* seen one which
1292 # would be indicated by None.
1293 #
1294 # For each argument there is the possibility that is_zero might give
1295 # True, False or None and likewise that is_infinite might give True,
1296 # False or None, giving 9 combinations. The True cases for is_zero and
1297 # is_infinite are mutually exclusive though so there are 3 main cases:
1298 #
1299 # - is_zero = True
1300 # - is_infinite = True
1301 # - is_zero and is_infinite are both either False or None
1302 #
1303 # At the end seen_zero and seen_infinite can be any of 9 combinations
1304 # of True/False/None. Unless one is False though we cannot return
1305 # anything except None:
1306 #
1307 # - is_zero=True needs seen_zero=True and seen_infinite=False
1308 # - is_zero=False needs seen_zero=False
1309 # - is_infinite=True needs seen_infinite=True and seen_zero=False
1310 # - is_infinite=False needs seen_infinite=False
1311 # - anything else gives both is_zero=None and is_infinite=None
1312 #
1313 # The loop only sets the flags to True or None and never back to False.
1314 # Hence as soon as neither flag is False we exit early returning None.
1315 # In particular as soon as we encounter a single arg that has
1316 # is_zero=is_infinite=None we exit. This is a common case since it is
1317 # the default assumptions for a Symbol and also the case for most
1318 # expressions containing such a symbol. The early exit gives a big
1319 # speedup for something like Mul(*symbols('x:1000')).is_zero.
1320 #
1321 seen_zero = seen_infinite = False
1323 for a in self.args:
1324 if a.is_zero:
1325 if seen_infinite is not False:
1326 return None, None
1327 seen_zero = True
1328 elif a.is_infinite:
1329 if seen_zero is not False:
1330 return None, None
1331 seen_infinite = True
1332 else:
1333 if seen_zero is False and a.is_zero is None:
1334 if seen_infinite is not False:
1335 return None, None
1336 seen_zero = None
1337 if seen_infinite is False and a.is_infinite is None:
1338 if seen_zero is not False:
1339 return None, None
1340 seen_infinite = None
1342 return seen_zero, seen_infinite
1344 def _eval_is_zero(self):
1345 # True iff any arg is zero and no arg is infinite but need to handle
1346 # three valued logic carefully.
1347 seen_zero, seen_infinite = self._eval_is_zero_infinite_helper()
1349 if seen_zero is False:
1350 return False
1351 elif seen_zero is True and seen_infinite is False:
1352 return True
1353 else:
1354 return None
1356 def _eval_is_infinite(self):
1357 # True iff any arg is infinite and no arg is zero but need to handle
1358 # three valued logic carefully.
1359 seen_zero, seen_infinite = self._eval_is_zero_infinite_helper()
1361 if seen_infinite is True and seen_zero is False:
1362 return True
1363 elif seen_infinite is False:
1364 return False
1365 else:
1366 return None
1368 # We do not need to implement _eval_is_finite because the assumptions
1369 # system can infer it from finite = not infinite.
1371 def _eval_is_rational(self):
1372 r = _fuzzy_group((a.is_rational for a in self.args), quick_exit=True)
1373 if r:
1374 return r
1375 elif r is False:
1376 # All args except one are rational
1377 if all(a.is_zero is False for a in self.args):
1378 return False
1380 def _eval_is_algebraic(self):
1381 r = _fuzzy_group((a.is_algebraic for a in self.args), quick_exit=True)
1382 if r:
1383 return r
1384 elif r is False:
1385 # All args except one are algebraic
1386 if all(a.is_zero is False for a in self.args):
1387 return False
1389 # without involving odd/even checks this code would suffice:
1390 #_eval_is_integer = lambda self: _fuzzy_group(
1391 # (a.is_integer for a in self.args), quick_exit=True)
1392 def _eval_is_integer(self):
1393 from sympy.ntheory.factor_ import trailing
1394 is_rational = self._eval_is_rational()
1395 if is_rational is False:
1396 return False
1398 numerators = []
1399 denominators = []
1400 unknown = False
1401 for a in self.args:
1402 hit = False
1403 if a.is_integer:
1404 if abs(a) is not S.One:
1405 numerators.append(a)
1406 elif a.is_Rational:
1407 n, d = a.as_numer_denom()
1408 if abs(n) is not S.One:
1409 numerators.append(n)
1410 if d is not S.One:
1411 denominators.append(d)
1412 elif a.is_Pow:
1413 b, e = a.as_base_exp()
1414 if not b.is_integer or not e.is_integer:
1415 hit = unknown = True
1416 if e.is_negative:
1417 denominators.append(2 if a is S.Half else
1418 Pow(a, S.NegativeOne))
1419 elif not hit:
1420 # int b and pos int e: a = b**e is integer
1421 assert not e.is_positive
1422 # for rational self and e equal to zero: a = b**e is 1
1423 assert not e.is_zero
1424 return # sign of e unknown -> self.is_integer unknown
1425 else:
1426 # x**2, 2**x, or x**y with x and y int-unknown -> unknown
1427 return
1428 else:
1429 return
1431 if not denominators and not unknown:
1432 return True
1434 allodd = lambda x: all(i.is_odd for i in x)
1435 alleven = lambda x: all(i.is_even for i in x)
1436 anyeven = lambda x: any(i.is_even for i in x)
1438 from .relational import is_gt
1439 if not numerators and denominators and all(
1440 is_gt(_, S.One) for _ in denominators):
1441 return False
1442 elif unknown:
1443 return
1444 elif allodd(numerators) and anyeven(denominators):
1445 return False
1446 elif anyeven(numerators) and denominators == [2]:
1447 return True
1448 elif alleven(numerators) and allodd(denominators
1449 ) and (Mul(*denominators, evaluate=False) - 1
1450 ).is_positive:
1451 return False
1452 if len(denominators) == 1:
1453 d = denominators[0]
1454 if d.is_Integer and d.is_even:
1455 # if minimal power of 2 in num vs den is not
1456 # negative then we have an integer
1457 if (Add(*[i.as_base_exp()[1] for i in
1458 numerators if i.is_even]) - trailing(d.p)
1459 ).is_nonnegative:
1460 return True
1461 if len(numerators) == 1:
1462 n = numerators[0]
1463 if n.is_Integer and n.is_even:
1464 # if minimal power of 2 in den vs num is positive
1465 # then we have have a non-integer
1466 if (Add(*[i.as_base_exp()[1] for i in
1467 denominators if i.is_even]) - trailing(n.p)
1468 ).is_positive:
1469 return False
1471 def _eval_is_polar(self):
1472 has_polar = any(arg.is_polar for arg in self.args)
1473 return has_polar and \
1474 all(arg.is_polar or arg.is_positive for arg in self.args)
1476 def _eval_is_extended_real(self):
1477 return self._eval_real_imag(True)
1479 def _eval_real_imag(self, real):
1480 zero = False
1481 t_not_re_im = None
1483 for t in self.args:
1484 if (t.is_complex or t.is_infinite) is False and t.is_extended_real is False:
1485 return False
1486 elif t.is_imaginary: # I
1487 real = not real
1488 elif t.is_extended_real: # 2
1489 if not zero:
1490 z = t.is_zero
1491 if not z and zero is False:
1492 zero = z
1493 elif z:
1494 if all(a.is_finite for a in self.args):
1495 return True
1496 return
1497 elif t.is_extended_real is False:
1498 # symbolic or literal like `2 + I` or symbolic imaginary
1499 if t_not_re_im:
1500 return # complex terms might cancel
1501 t_not_re_im = t
1502 elif t.is_imaginary is False: # symbolic like `2` or `2 + I`
1503 if t_not_re_im:
1504 return # complex terms might cancel
1505 t_not_re_im = t
1506 else:
1507 return
1509 if t_not_re_im:
1510 if t_not_re_im.is_extended_real is False:
1511 if real: # like 3
1512 return zero # 3*(smthng like 2 + I or i) is not real
1513 if t_not_re_im.is_imaginary is False: # symbolic 2 or 2 + I
1514 if not real: # like I
1515 return zero # I*(smthng like 2 or 2 + I) is not real
1516 elif zero is False:
1517 return real # can't be trumped by 0
1518 elif real:
1519 return real # doesn't matter what zero is
1521 def _eval_is_imaginary(self):
1522 if all(a.is_zero is False and a.is_finite for a in self.args):
1523 return self._eval_real_imag(False)
1525 def _eval_is_hermitian(self):
1526 return self._eval_herm_antiherm(True)
1528 def _eval_is_antihermitian(self):
1529 return self._eval_herm_antiherm(False)
1531 def _eval_herm_antiherm(self, herm):
1532 for t in self.args:
1533 if t.is_hermitian is None or t.is_antihermitian is None:
1534 return
1535 if t.is_hermitian:
1536 continue
1537 elif t.is_antihermitian:
1538 herm = not herm
1539 else:
1540 return
1542 if herm is not False:
1543 return herm
1545 is_zero = self._eval_is_zero()
1546 if is_zero:
1547 return True
1548 elif is_zero is False:
1549 return herm
1551 def _eval_is_irrational(self):
1552 for t in self.args:
1553 a = t.is_irrational
1554 if a:
1555 others = list(self.args)
1556 others.remove(t)
1557 if all((x.is_rational and fuzzy_not(x.is_zero)) is True for x in others):
1558 return True
1559 return
1560 if a is None:
1561 return
1562 if all(x.is_real for x in self.args):
1563 return False
1565 def _eval_is_extended_positive(self):
1566 """Return True if self is positive, False if not, and None if it
1567 cannot be determined.
1569 Explanation
1570 ===========
1572 This algorithm is non-recursive and works by keeping track of the
1573 sign which changes when a negative or nonpositive is encountered.
1574 Whether a nonpositive or nonnegative is seen is also tracked since
1575 the presence of these makes it impossible to return True, but
1576 possible to return False if the end result is nonpositive. e.g.
1578 pos * neg * nonpositive -> pos or zero -> None is returned
1579 pos * neg * nonnegative -> neg or zero -> False is returned
1580 """
1581 return self._eval_pos_neg(1)
1583 def _eval_pos_neg(self, sign):
1584 saw_NON = saw_NOT = False
1585 for t in self.args:
1586 if t.is_extended_positive:
1587 continue
1588 elif t.is_extended_negative:
1589 sign = -sign
1590 elif t.is_zero:
1591 if all(a.is_finite for a in self.args):
1592 return False
1593 return
1594 elif t.is_extended_nonpositive:
1595 sign = -sign
1596 saw_NON = True
1597 elif t.is_extended_nonnegative:
1598 saw_NON = True
1599 # FIXME: is_positive/is_negative is False doesn't take account of
1600 # Symbol('x', infinite=True, extended_real=True) which has
1601 # e.g. is_positive is False but has uncertain sign.
1602 elif t.is_positive is False:
1603 sign = -sign
1604 if saw_NOT:
1605 return
1606 saw_NOT = True
1607 elif t.is_negative is False:
1608 if saw_NOT:
1609 return
1610 saw_NOT = True
1611 else:
1612 return
1613 if sign == 1 and saw_NON is False and saw_NOT is False:
1614 return True
1615 if sign < 0:
1616 return False
1618 def _eval_is_extended_negative(self):
1619 return self._eval_pos_neg(-1)
1621 def _eval_is_odd(self):
1622 is_integer = self._eval_is_integer()
1623 if is_integer is not True:
1624 return is_integer
1626 from sympy.simplify.radsimp import fraction
1627 n, d = fraction(self)
1628 if d.is_Integer and d.is_even:
1629 from sympy.ntheory.factor_ import trailing
1630 # if minimal power of 2 in num vs den is
1631 # positive then we have an even number
1632 if (Add(*[i.as_base_exp()[1] for i in
1633 Mul.make_args(n) if i.is_even]) - trailing(d.p)
1634 ).is_positive:
1635 return False
1636 return
1637 r, acc = True, 1
1638 for t in self.args:
1639 if abs(t) is S.One:
1640 continue
1641 if t.is_even:
1642 return False
1643 if r is False:
1644 pass
1645 elif acc != 1 and (acc + t).is_odd:
1646 r = False
1647 elif t.is_even is None:
1648 r = None
1649 acc = t
1650 return r
1652 def _eval_is_even(self):
1653 from sympy.simplify.radsimp import fraction
1654 n, d = fraction(self)
1655 if n.is_Integer and n.is_even:
1656 # if minimal power of 2 in den vs num is not
1657 # negative then this is not an integer and
1658 # can't be even
1659 from sympy.ntheory.factor_ import trailing
1660 if (Add(*[i.as_base_exp()[1] for i in
1661 Mul.make_args(d) if i.is_even]) - trailing(n.p)
1662 ).is_nonnegative:
1663 return False
1665 def _eval_is_composite(self):
1666 """
1667 Here we count the number of arguments that have a minimum value
1668 greater than two.
1669 If there are more than one of such a symbol then the result is composite.
1670 Else, the result cannot be determined.
1671 """
1672 number_of_args = 0 # count of symbols with minimum value greater than one
1673 for arg in self.args:
1674 if not (arg.is_integer and arg.is_positive):
1675 return None
1676 if (arg-1).is_positive:
1677 number_of_args += 1
1679 if number_of_args > 1:
1680 return True
1682 def _eval_subs(self, old, new):
1683 from sympy.functions.elementary.complexes import sign
1684 from sympy.ntheory.factor_ import multiplicity
1685 from sympy.simplify.powsimp import powdenest
1686 from sympy.simplify.radsimp import fraction
1688 if not old.is_Mul:
1689 return None
1691 # try keep replacement literal so -2*x doesn't replace 4*x
1692 if old.args[0].is_Number and old.args[0] < 0:
1693 if self.args[0].is_Number:
1694 if self.args[0] < 0:
1695 return self._subs(-old, -new)
1696 return None
1698 def base_exp(a):
1699 # if I and -1 are in a Mul, they get both end up with
1700 # a -1 base (see issue 6421); all we want here are the
1701 # true Pow or exp separated into base and exponent
1702 from sympy.functions.elementary.exponential import exp
1703 if a.is_Pow or isinstance(a, exp):
1704 return a.as_base_exp()
1705 return a, S.One
1707 def breakup(eq):
1708 """break up powers of eq when treated as a Mul:
1709 b**(Rational*e) -> b**e, Rational
1710 commutatives come back as a dictionary {b**e: Rational}
1711 noncommutatives come back as a list [(b**e, Rational)]
1712 """
1714 (c, nc) = (defaultdict(int), [])
1715 for a in Mul.make_args(eq):
1716 a = powdenest(a)
1717 (b, e) = base_exp(a)
1718 if e is not S.One:
1719 (co, _) = e.as_coeff_mul()
1720 b = Pow(b, e/co)
1721 e = co
1722 if a.is_commutative:
1723 c[b] += e
1724 else:
1725 nc.append([b, e])
1726 return (c, nc)
1728 def rejoin(b, co):
1729 """
1730 Put rational back with exponent; in general this is not ok, but
1731 since we took it from the exponent for analysis, it's ok to put
1732 it back.
1733 """
1735 (b, e) = base_exp(b)
1736 return Pow(b, e*co)
1738 def ndiv(a, b):
1739 """if b divides a in an extractive way (like 1/4 divides 1/2
1740 but not vice versa, and 2/5 does not divide 1/3) then return
1741 the integer number of times it divides, else return 0.
1742 """
1743 if not b.q % a.q or not a.q % b.q:
1744 return int(a/b)
1745 return 0
1747 # give Muls in the denominator a chance to be changed (see issue 5651)
1748 # rv will be the default return value
1749 rv = None
1750 n, d = fraction(self)
1751 self2 = self
1752 if d is not S.One:
1753 self2 = n._subs(old, new)/d._subs(old, new)
1754 if not self2.is_Mul:
1755 return self2._subs(old, new)
1756 if self2 != self:
1757 rv = self2
1759 # Now continue with regular substitution.
1761 # handle the leading coefficient and use it to decide if anything
1762 # should even be started; we always know where to find the Rational
1763 # so it's a quick test
1765 co_self = self2.args[0]
1766 co_old = old.args[0]
1767 co_xmul = None
1768 if co_old.is_Rational and co_self.is_Rational:
1769 # if coeffs are the same there will be no updating to do
1770 # below after breakup() step; so skip (and keep co_xmul=None)
1771 if co_old != co_self:
1772 co_xmul = co_self.extract_multiplicatively(co_old)
1773 elif co_old.is_Rational:
1774 return rv
1776 # break self and old into factors
1778 (c, nc) = breakup(self2)
1779 (old_c, old_nc) = breakup(old)
1781 # update the coefficients if we had an extraction
1782 # e.g. if co_self were 2*(3/35*x)**2 and co_old = 3/5
1783 # then co_self in c is replaced by (3/5)**2 and co_residual
1784 # is 2*(1/7)**2
1786 if co_xmul and co_xmul.is_Rational and abs(co_old) != 1:
1787 mult = S(multiplicity(abs(co_old), co_self))
1788 c.pop(co_self)
1789 if co_old in c:
1790 c[co_old] += mult
1791 else:
1792 c[co_old] = mult
1793 co_residual = co_self/co_old**mult
1794 else:
1795 co_residual = 1
1797 # do quick tests to see if we can't succeed
1799 ok = True
1800 if len(old_nc) > len(nc):
1801 # more non-commutative terms
1802 ok = False
1803 elif len(old_c) > len(c):
1804 # more commutative terms
1805 ok = False
1806 elif {i[0] for i in old_nc}.difference({i[0] for i in nc}):
1807 # unmatched non-commutative bases
1808 ok = False
1809 elif set(old_c).difference(set(c)):
1810 # unmatched commutative terms
1811 ok = False
1812 elif any(sign(c[b]) != sign(old_c[b]) for b in old_c):
1813 # differences in sign
1814 ok = False
1815 if not ok:
1816 return rv
1818 if not old_c:
1819 cdid = None
1820 else:
1821 rat = []
1822 for (b, old_e) in old_c.items():
1823 c_e = c[b]
1824 rat.append(ndiv(c_e, old_e))
1825 if not rat[-1]:
1826 return rv
1827 cdid = min(rat)
1829 if not old_nc:
1830 ncdid = None
1831 for i in range(len(nc)):
1832 nc[i] = rejoin(*nc[i])
1833 else:
1834 ncdid = 0 # number of nc replacements we did
1835 take = len(old_nc) # how much to look at each time
1836 limit = cdid or S.Infinity # max number that we can take
1837 failed = [] # failed terms will need subs if other terms pass
1838 i = 0
1839 while limit and i + take <= len(nc):
1840 hit = False
1842 # the bases must be equivalent in succession, and
1843 # the powers must be extractively compatible on the
1844 # first and last factor but equal in between.
1846 rat = []
1847 for j in range(take):
1848 if nc[i + j][0] != old_nc[j][0]:
1849 break
1850 elif j == 0:
1851 rat.append(ndiv(nc[i + j][1], old_nc[j][1]))
1852 elif j == take - 1:
1853 rat.append(ndiv(nc[i + j][1], old_nc[j][1]))
1854 elif nc[i + j][1] != old_nc[j][1]:
1855 break
1856 else:
1857 rat.append(1)
1858 j += 1
1859 else:
1860 ndo = min(rat)
1861 if ndo:
1862 if take == 1:
1863 if cdid:
1864 ndo = min(cdid, ndo)
1865 nc[i] = Pow(new, ndo)*rejoin(nc[i][0],
1866 nc[i][1] - ndo*old_nc[0][1])
1867 else:
1868 ndo = 1
1870 # the left residual
1872 l = rejoin(nc[i][0], nc[i][1] - ndo*
1873 old_nc[0][1])
1875 # eliminate all middle terms
1877 mid = new
1879 # the right residual (which may be the same as the middle if take == 2)
1881 ir = i + take - 1
1882 r = (nc[ir][0], nc[ir][1] - ndo*
1883 old_nc[-1][1])
1884 if r[1]:
1885 if i + take < len(nc):
1886 nc[i:i + take] = [l*mid, r]
1887 else:
1888 r = rejoin(*r)
1889 nc[i:i + take] = [l*mid*r]
1890 else:
1892 # there was nothing left on the right
1894 nc[i:i + take] = [l*mid]
1896 limit -= ndo
1897 ncdid += ndo
1898 hit = True
1899 if not hit:
1901 # do the subs on this failing factor
1903 failed.append(i)
1904 i += 1
1905 else:
1907 if not ncdid:
1908 return rv
1910 # although we didn't fail, certain nc terms may have
1911 # failed so we rebuild them after attempting a partial
1912 # subs on them
1914 failed.extend(range(i, len(nc)))
1915 for i in failed:
1916 nc[i] = rejoin(*nc[i]).subs(old, new)
1918 # rebuild the expression
1920 if cdid is None:
1921 do = ncdid
1922 elif ncdid is None:
1923 do = cdid
1924 else:
1925 do = min(ncdid, cdid)
1927 margs = []
1928 for b in c:
1929 if b in old_c:
1931 # calculate the new exponent
1933 e = c[b] - old_c[b]*do
1934 margs.append(rejoin(b, e))
1935 else:
1936 margs.append(rejoin(b.subs(old, new), c[b]))
1937 if cdid and not ncdid:
1939 # in case we are replacing commutative with non-commutative,
1940 # we want the new term to come at the front just like the
1941 # rest of this routine
1943 margs = [Pow(new, cdid)] + margs
1944 return co_residual*self2.func(*margs)*self2.func(*nc)
1946 def _eval_nseries(self, x, n, logx, cdir=0):
1947 from .function import PoleError
1948 from sympy.functions.elementary.integers import ceiling
1949 from sympy.series.order import Order
1951 def coeff_exp(term, x):
1952 lt = term.as_coeff_exponent(x)
1953 if lt[0].has(x):
1954 try:
1955 lt = term.leadterm(x)
1956 except ValueError:
1957 return term, S.Zero
1958 return lt
1960 ords = []
1962 try:
1963 for t in self.args:
1964 coeff, exp = t.leadterm(x)
1965 if not coeff.has(x):
1966 ords.append((t, exp))
1967 else:
1968 raise ValueError
1970 n0 = sum(t[1] for t in ords if t[1].is_number)
1971 facs = []
1972 for t, m in ords:
1973 n1 = ceiling(n - n0 + (m if m.is_number else 0))
1974 s = t.nseries(x, n=n1, logx=logx, cdir=cdir)
1975 ns = s.getn()
1976 if ns is not None:
1977 if ns < n1: # less than expected
1978 n -= n1 - ns # reduce n
1979 facs.append(s)
1981 except (ValueError, NotImplementedError, TypeError, AttributeError, PoleError):
1982 n0 = sympify(sum(t[1] for t in ords if t[1].is_number))
1983 if n0.is_nonnegative:
1984 n0 = S.Zero
1985 facs = [t.nseries(x, n=ceiling(n-n0), logx=logx, cdir=cdir) for t in self.args]
1986 from sympy.simplify.powsimp import powsimp
1987 res = powsimp(self.func(*facs).expand(), combine='exp', deep=True)
1988 if res.has(Order):
1989 res += Order(x**n, x)
1990 return res
1992 res = S.Zero
1993 ords2 = [Add.make_args(factor) for factor in facs]
1995 for fac in product(*ords2):
1996 ords3 = [coeff_exp(term, x) for term in fac]
1997 coeffs, powers = zip(*ords3)
1998 power = sum(powers)
1999 if (power - n).is_negative:
2000 res += Mul(*coeffs)*(x**power)
2002 def max_degree(e, x):
2003 if e is x:
2004 return S.One
2005 if e.is_Atom:
2006 return S.Zero
2007 if e.is_Add:
2008 return max(max_degree(a, x) for a in e.args)
2009 if e.is_Mul:
2010 return Add(*[max_degree(a, x) for a in e.args])
2011 if e.is_Pow:
2012 return max_degree(e.base, x)*e.exp
2013 return S.Zero
2015 if self.is_polynomial(x):
2016 from sympy.polys.polyerrors import PolynomialError
2017 from sympy.polys.polytools import degree
2018 try:
2019 if max_degree(self, x) >= n or degree(self, x) != degree(res, x):
2020 res += Order(x**n, x)
2021 except PolynomialError:
2022 pass
2023 else:
2024 return res
2026 if res != self:
2027 if (self - res).subs(x, 0) == S.Zero and n > 0:
2028 lt = self._eval_as_leading_term(x, logx=logx, cdir=cdir)
2029 if lt == S.Zero:
2030 return res
2031 res += Order(x**n, x)
2032 return res
2034 def _eval_as_leading_term(self, x, logx=None, cdir=0):
2035 return self.func(*[t.as_leading_term(x, logx=logx, cdir=cdir) for t in self.args])
2037 def _eval_conjugate(self):
2038 return self.func(*[t.conjugate() for t in self.args])
2040 def _eval_transpose(self):
2041 return self.func(*[t.transpose() for t in self.args[::-1]])
2043 def _eval_adjoint(self):
2044 return self.func(*[t.adjoint() for t in self.args[::-1]])
2046 def as_content_primitive(self, radical=False, clear=True):
2047 """Return the tuple (R, self/R) where R is the positive Rational
2048 extracted from self.
2050 Examples
2051 ========
2053 >>> from sympy import sqrt
2054 >>> (-3*sqrt(2)*(2 - 2*sqrt(2))).as_content_primitive()
2055 (6, -sqrt(2)*(1 - sqrt(2)))
2057 See docstring of Expr.as_content_primitive for more examples.
2058 """
2060 coef = S.One
2061 args = []
2062 for a in self.args:
2063 c, p = a.as_content_primitive(radical=radical, clear=clear)
2064 coef *= c
2065 if p is not S.One:
2066 args.append(p)
2067 # don't use self._from_args here to reconstruct args
2068 # since there may be identical args now that should be combined
2069 # e.g. (2+2*x)*(3+3*x) should be (6, (1 + x)**2) not (6, (1+x)*(1+x))
2070 return coef, self.func(*args)
2072 def as_ordered_factors(self, order=None):
2073 """Transform an expression into an ordered list of factors.
2075 Examples
2076 ========
2078 >>> from sympy import sin, cos
2079 >>> from sympy.abc import x, y
2081 >>> (2*x*y*sin(x)*cos(x)).as_ordered_factors()
2082 [2, x, y, sin(x), cos(x)]
2084 """
2085 cpart, ncpart = self.args_cnc()
2086 cpart.sort(key=lambda expr: expr.sort_key(order=order))
2087 return cpart + ncpart
2089 @property
2090 def _sorted_args(self):
2091 return tuple(self.as_ordered_factors())
2093mul = AssocOpDispatcher('mul')
2096def prod(a, start=1):
2097 """Return product of elements of a. Start with int 1 so if only
2098 ints are included then an int result is returned.
2100 Examples
2101 ========
2103 >>> from sympy import prod, S
2104 >>> prod(range(3))
2105 0
2106 >>> type(_) is int
2107 True
2108 >>> prod([S(2), 3])
2109 6
2110 >>> _.is_Integer
2111 True
2113 You can start the product at something other than 1:
2115 >>> prod([1, 2], 3)
2116 6
2118 """
2119 return reduce(operator.mul, a, start)
2122def _keep_coeff(coeff, factors, clear=True, sign=False):
2123 """Return ``coeff*factors`` unevaluated if necessary.
2125 If ``clear`` is False, do not keep the coefficient as a factor
2126 if it can be distributed on a single factor such that one or
2127 more terms will still have integer coefficients.
2129 If ``sign`` is True, allow a coefficient of -1 to remain factored out.
2131 Examples
2132 ========
2134 >>> from sympy.core.mul import _keep_coeff
2135 >>> from sympy.abc import x, y
2136 >>> from sympy import S
2138 >>> _keep_coeff(S.Half, x + 2)
2139 (x + 2)/2
2140 >>> _keep_coeff(S.Half, x + 2, clear=False)
2141 x/2 + 1
2142 >>> _keep_coeff(S.Half, (x + 2)*y, clear=False)
2143 y*(x + 2)/2
2144 >>> _keep_coeff(S(-1), x + y)
2145 -x - y
2146 >>> _keep_coeff(S(-1), x + y, sign=True)
2147 -(x + y)
2148 """
2149 if not coeff.is_Number:
2150 if factors.is_Number:
2151 factors, coeff = coeff, factors
2152 else:
2153 return coeff*factors
2154 if factors is S.One:
2155 return coeff
2156 if coeff is S.One:
2157 return factors
2158 elif coeff is S.NegativeOne and not sign:
2159 return -factors
2160 elif factors.is_Add:
2161 if not clear and coeff.is_Rational and coeff.q != 1:
2162 args = [i.as_coeff_Mul() for i in factors.args]
2163 args = [(_keep_coeff(c, coeff), m) for c, m in args]
2164 if any(c.is_Integer for c, _ in args):
2165 return Add._from_args([Mul._from_args(
2166 i[1:] if i[0] == 1 else i) for i in args])
2167 return Mul(coeff, factors, evaluate=False)
2168 elif factors.is_Mul:
2169 margs = list(factors.args)
2170 if margs[0].is_Number:
2171 margs[0] *= coeff
2172 if margs[0] == 1:
2173 margs.pop(0)
2174 else:
2175 margs.insert(0, coeff)
2176 return Mul._from_args(margs)
2177 else:
2178 m = coeff*factors
2179 if m.is_Number and not factors.is_Number:
2180 m = Mul._from_args((coeff, factors))
2181 return m
2183def expand_2arg(e):
2184 def do(e):
2185 if e.is_Mul:
2186 c, r = e.as_coeff_Mul()
2187 if c.is_Number and r.is_Add:
2188 return _unevaluated_Add(*[c*ri for ri in r.args])
2189 return e
2190 return bottom_up(e, do)
2193from .numbers import Rational
2194from .power import Pow
2195from .add import Add, _unevaluated_Add