Hot-keys on this page
r m x p toggle line displays
j k next/prev highlighted chunk
0 (zero) top of page
1 (one) first highlighted chunk
from __future__ import print_function, division
from collections import defaultdict
from sympy.core.function import expand_log, count_ops from sympy.core import sympify, Basic, Dummy, S, Add, Mul, Pow, expand_mul, factor_terms from sympy.core.compatibility import ordered, default_sort_key, reduce from sympy.core.numbers import Integer, Rational from sympy.core.mul import prod, _keep_coeff from sympy.core.rules import Transform from sympy.functions import exp_polar, exp, log, root, polarify, unpolarify from sympy.polys import lcm, gcd from sympy.ntheory.factor_ import multiplicity
def powsimp(expr, deep=False, combine='all', force=False, measure=count_ops): """ reduces expression by combining powers with similar bases and exponents.
Notes =====
If deep is True then powsimp() will also simplify arguments of functions. By default deep is set to False.
If force is True then bases will be combined without checking for assumptions, e.g. sqrt(x)*sqrt(y) -> sqrt(x*y) which is not true if x and y are both negative.
You can make powsimp() only combine bases or only combine exponents by changing combine='base' or combine='exp'. By default, combine='all', which does both. combine='base' will only combine::
a a a 2x x x * y => (x*y) as well as things like 2 => 4
and combine='exp' will only combine ::
a b (a + b) x * x => x
combine='exp' will strictly only combine exponents in the way that used to be automatic. Also use deep=True if you need the old behavior.
When combine='all', 'exp' is evaluated first. Consider the first example below for when there could be an ambiguity relating to this. This is done so things like the second example can be completely combined. If you want 'base' combined first, do something like powsimp(powsimp(expr, combine='base'), combine='exp').
Examples ========
>>> from sympy import powsimp, exp, log, symbols >>> from sympy.abc import x, y, z, n >>> powsimp(x**y*x**z*y**z, combine='all') x**(y + z)*y**z >>> powsimp(x**y*x**z*y**z, combine='exp') x**(y + z)*y**z >>> powsimp(x**y*x**z*y**z, combine='base', force=True) x**y*(x*y)**z
>>> powsimp(x**z*x**y*n**z*n**y, combine='all', force=True) (n*x)**(y + z) >>> powsimp(x**z*x**y*n**z*n**y, combine='exp') n**(y + z)*x**(y + z) >>> powsimp(x**z*x**y*n**z*n**y, combine='base', force=True) (n*x)**y*(n*x)**z
>>> x, y = symbols('x y', positive=True) >>> powsimp(log(exp(x)*exp(y))) log(exp(x)*exp(y)) >>> powsimp(log(exp(x)*exp(y)), deep=True) x + y
Radicals with Mul bases will be combined if combine='exp'
>>> from sympy import sqrt, Mul >>> x, y = symbols('x y')
Two radicals are automatically joined through Mul:
>>> a=sqrt(x*sqrt(y)) >>> a*a**3 == a**4 True
But if an integer power of that radical has been autoexpanded then Mul does not join the resulting factors:
>>> a**4 # auto expands to a Mul, no longer a Pow x**2*y >>> _*a # so Mul doesn't combine them x**2*y*sqrt(x*sqrt(y)) >>> powsimp(_) # but powsimp will (x*sqrt(y))**(5/2) >>> powsimp(x*y*a) # but won't when doing so would violate assumptions x*y*sqrt(x*sqrt(y))
"""
expr.is_Atom or expr in (exp_polar(0), exp_polar(1)))):
# handle the Mul # Collect base/exp data, while maintaining order in the # non-commutative parts of the product # don't let smthg like sqrt(x**a) split into x**a, 1/2 # or else it will be joined as x**(a/2) later b, e = b**e, S.One else: # This is the logic that combines exponents for equal, # but non-commutative bases: A**x*A**y == A**(x+y). if nc_part: b1, e1 = nc_part[-1].as_base_exp() b2, e2 = term.as_base_exp() if (b1 == b2 and e1.is_commutative and e2.is_commutative): nc_part[-1] = Pow(b1, Add(e1, e2)) continue nc_part.append(term)
# add up exponents of common bases # allow 2**x/4 -> 2**(x - 2); don't do this when b and e are # Numbers since autoevaluation will undo it, e.g. # 2**(1/3)/4 -> 2**(1/3 - 2) -> 2**(1/3)/4 coeff is not S.One and b not in (S.One, S.NegativeOne)): m = multiplicity(abs(b), abs(coeff)) if m: e.append(m) coeff /= b**m else:
# convert to plain dictionary
# check for base and inverted base pairs continue else: skip.add(binv) e = c_powers.pop(binv) c_powers[b] -= e
# check for base and negated base pairs if (b.is_positive in (0, 1) or e.is_integer): c_powers[-b] += c_powers.pop(b) if _n in c_powers: c_powers[_n] += e else: c_powers[_n] = e
# filter c_powers and convert to a list
# ============================================================== # check for Mul bases of Rational powers that can be combined with # separated bases, e.g. x*sqrt(x*y)*sqrt(x*sqrt(x*y)) -> # (x*sqrt(x*y))**(3/2) # ---------------- helper functions
'''Return Rational part of x's exponent as it appears in the bkey. '''
'''Return (b**s, c.q), c.p where e -> c*s. If e is not given then it will be taken by using as_base_exp() on the input b. e.g. x**3/2 -> (x, 2), 3 x**y -> (x**y, 1), 1 x**(2*y/3) -> (x**y, 3), 2 exp(x/2) -> (exp(a), 2), 1
''' else: return (b, Integer(c.q)), m*Integer(c.p) else: else:
'''Decide what to do with base, b. If its exponent is now an integer multiple of the Rational denominator, then remove it and put the factors of its base in the common_b dictionary or update the existing bases if necessary. If it has been zeroed out, simply remove the base. ''' newe, r = divmod(common_b[b], b[1]) if not r: common_b.pop(b) if newe: for m in Mul.make_args(b[0]**newe): b, e = bkey(m) if b not in common_b: common_b[b] = 0 common_b[b] += e if b[1] != 1: bases.append(b) # ---------------- end of helper functions
# assemble a dictionary of the factors having a Rational power common_b[b] = common_b[b] + e else: continue
# find the number of extractions possible # e.g. [(1, 2), (2, 2)] -> min(2/1, 2/2) -> 1 min1 = ee[0][1]/ee[0][0] for i in range(len(ee)): rat = ee[i][1]/ee[i][0] if rat < 1: break min1 = min(min1, rat) else: # update base factor counts # e.g. if ee = [(2, 5), (3, 6)] then min1 = 2 # and the new base counts will be 5-2*2 and 6-2*3 for i in range(len(bb)): common_b[bb[i]] -= min1*ee[i][0] update(bb[i]) # update the count of the base # e.g. x**2*y*sqrt(x*sqrt(y)) the count of x*sqrt(y) # will increase by 4 to give bkey (x*sqrt(y), 2, 5) common_b[base] += min1*qstart*exponent or len(common_b) == 1 # nothing left to join with or all(k[1] == 1 for k in common_b) # no rad's in common_b ): # see what we can exponentiate base by to remove any radicals # so we know what to search for # e.g. if base were x**(1/2)*y**(1/3) then we should # exponentiate by 6 and look for powers of x and y in the ratio # of 2 to 3 # this base no longer can find anything to join with and # since it was longer than any other we are done with it
# update c_powers and get ready to continue with powsimp # there may be terms still in common_b that were bases that were # identified as needing processing, so remove those, too q is not S.One and not b.exp.is_Rational: b, be = b.as_base_exp() b = b**(be/q) else: # ==============================================================
# rebuild the expression else: recurse(newexpr, combine='base')
# Build c_powers and nc_part. These must both be lists not # dicts because exp's are not combined. else: # This is the logic that combines bases that are # different and non-commutative, but with equal and # commutative exponents: A**x*B**x == (A*B)**x. if nc_part: b1, e1 = nc_part[-1].as_base_exp() b2, e2 = term.as_base_exp() if (e1 == e2 and e2.is_commutative): nc_part[-1] = Pow(b1*b2, e1) continue nc_part.append(term)
# Pull out numerical coefficients from exponent if assumptions allow # e.g., 2**(2*x) => 4**x
# Combine bases whenever they have the same exponent and # assumptions allow # first gather the potential bases under the common exponent e = recurse(e)
# Merge back in the results of the above to form a new product
# calculate the new base for e
else: # see which ones can be joined nonneg.append( bi) # polar can be treated like non-negative else: # a single neg or a single unk can join the rest # their negative signs cancel in groups of 2*q if we know # that e = p/q else we have to treat them as unknown israt = False if e.is_Rational: israt = True else: p, d = e.as_numer_denom() if p.is_integer and d.is_integer: israt = True if israt: neg = [-w for w in neg] unk.extend([S.NegativeOne]*len(neg)) else: unk.extend(neg) neg = [] del israt
# these shouldn't be joined # here is a new joined base # if there are positive parts they will just get separated # again unless some change is made
# return the number of terms of this expression # when multiplied out -- assuming no joining of terms new_base = factor_terms(xnew_base)
# break out the powers from c_powers now
# we're done
else: raise ValueError("combine must be one of ('all', 'exp', 'base').")
def powdenest(eq, force=False, polar=False): r""" Collect exponents on powers as assumptions allow.
Given ``(bb**be)**e``, this can be simplified as follows: * if ``bb`` is positive, or * ``e`` is an integer, or * ``|be| < 1`` then this simplifies to ``bb**(be*e)``
Given a product of powers raised to a power, ``(bb1**be1 * bb2**be2...)**e``, simplification can be done as follows:
- if e is positive, the gcd of all bei can be joined with e; - all non-negative bb can be separated from those that are negative and their gcd can be joined with e; autosimplification already handles this separation. - integer factors from powers that have integers in the denominator of the exponent can be removed from any term and the gcd of such integers can be joined with e
Setting ``force`` to True will make symbols that are not explicitly negative behave as though they are positive, resulting in more denesting.
Setting ``polar`` to True will do simplifications on the Riemann surface of the logarithm, also resulting in more denestings.
When there are sums of logs in exp() then a product of powers may be obtained e.g. ``exp(3*(log(a) + 2*log(b)))`` - > ``a**3*b**6``.
Examples ========
>>> from sympy.abc import a, b, x, y, z >>> from sympy import Symbol, exp, log, sqrt, symbols, powdenest
>>> powdenest((x**(2*a/3))**(3*x)) (x**(2*a/3))**(3*x) >>> powdenest(exp(3*x*log(2))) 2**(3*x)
Assumptions may prevent expansion:
>>> powdenest(sqrt(x**2)) sqrt(x**2)
>>> p = symbols('p', positive=True) >>> powdenest(sqrt(p**2)) p
No other expansion is done.
>>> i, j = symbols('i,j', integer=True) >>> powdenest((x**x)**(i + j)) # -X-> (x**x)**i*(x**x)**j x**(x*(i + j))
But exp() will be denested by moving all non-log terms outside of the function; this may result in the collapsing of the exp to a power with a different base:
>>> powdenest(exp(3*y*log(x))) x**(3*y) >>> powdenest(exp(y*(log(a) + log(b)))) (a*b)**y >>> powdenest(exp(3*(log(a) + log(b)))) a**3*b**3
If assumptions allow, symbols can also be moved to the outermost exponent:
>>> i = Symbol('i', integer=True) >>> powdenest(((x**(2*i))**(3*y))**x) ((x**(2*i))**(3*y))**x >>> powdenest(((x**(2*i))**(3*y))**x, force=True) x**(6*i*x*y)
>>> powdenest(((x**(2*a/3))**(3*y/i))**x) ((x**(2*a/3))**(3*y/i))**x >>> powdenest((x**(2*i)*y**(4*i))**z, force=True) (x*y**2)**(2*i*z)
>>> n = Symbol('n', negative=True)
>>> powdenest((x**i)**y, force=True) x**(i*y) >>> powdenest((n**i)**x, force=True) (n**i)**x
"""
eq, rep = posify(eq) return powdenest(eq, force=False).xreplace(rep)
eq, rep = polarify(eq) return unpolarify(powdenest(unpolarify(eq, exponents_only=True)), rep)
_denest_pow, filter=lambda m: m.is_Pow or m.func is exp))
_y = Dummy('y')
def _denest_pow(eq): """ Denest powers.
This is a helper function for powdenest that performs the actual transformation. """
new = b._eval_power(e) if new is not None: eq = new b, e = new.as_base_exp()
# denest exp with log terms in exponent logs.append(ei) else:
b.is_Rational and b.q != 1 or b.is_positive):
# denest eq which is either pos**e or Pow**e or Mul**e or # Mul(b1**e1, b2**e2)
# handle polar numbers specially polars.append(bb.as_base_exp()) else: return Pow(polars[0][0], polars[0][1]*e)*powdenest(Mul(*nonpolars)**e) return Mul(*[powdenest(bb**(ee*e)) for (bb, ee) in polars]) \ *powdenest(Mul(*nonpolars)**e)
# use log to see if there is a power here c, logb = logb.args e *= c base = logb.args[0] return Pow(base, e)
# if b is not a Mul or any factor is an atom then there is nothing to do
# let log handle the case of the base of the argument being a Mul, e.g. # sqrt(x**(2*i)*y**(6*i)) -> x**i*y**(3**i) if x and y are positive; we # will take the log, expand it, and then factor out the common powers that # now appear as coefficient. We do this manually since terms_gcd pulls out # fractions, terms_gcd(x+x*y/2) -> x*(y + 2)/2 and we don't want the 1/2; # gcd won't pull out numerators from a fraction: gcd(3*x, 9*x/2) -> x but # we want 3*x. Neither work with noncommutatives.
def nc_gcd(aa, bb): a, b = [i.as_coeff_Mul() for i in [aa, bb]] c = gcd(a[0], b[0]).as_numer_denom()[0] g = Mul(*(a[1].args_cnc(cset=True)[0] & b[1].args_cnc(cset=True)[0])) return _keep_coeff(c, g)
glogb = expand_log(log(b)) if glogb.is_Add: args = glogb.args g = reduce(nc_gcd, args) if g != 1: cg, rg = g.as_coeff_Mul() glogb = _keep_coeff(cg, rg*Add(*[a/g for a in args]))
# now put the log back together again if glogb.func is log or not glogb.is_Mul: if glogb.args[0].is_Pow or glogb.args[0].func is exp: glogb = _denest_pow(glogb.args[0]) if (abs(glogb.exp) < 1) == True: return Pow(glogb.base, glogb.exp*e) return eq
# the log(b) was a Mul so join any adds with logcombine add = [] other = [] for a in glogb.args: if a.is_Add: add.append(a) else: other.append(a) return Pow(exp(logcombine(Mul(*add))), e*Mul(*other)) |