本文整理汇总了Python中sympy.core.mul.Mul类的典型用法代码示例。如果您正苦于以下问题:Python Mul类的具体用法?Python Mul怎么用?Python Mul使用的例子?那么恭喜您, 这里精选的类代码示例或许可以为您提供帮助。
在下文中一共展示了Mul类的15个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的Python代码示例。
示例1: factor_terms
def factor_terms(expr, radical=False):
"""Remove common factors from terms in all arguments without
changing the underlying structure of the expr. No expansion or
simplification (and no processing of non-commutatives) is performed.
If radical=True then a radical common to all terms will be factored
out of any Add sub-expressions of the expr.
Examples
========
>>> from sympy import factor_terms, Symbol
>>> from sympy.abc import x, y
>>> factor_terms(x + x*(2 + 4*y)**3)
x*(8*(2*y + 1)**3 + 1)
>>> A = Symbol('A', commutative=False)
>>> factor_terms(x*A + x*A + x*y*A)
x*(y*A + 2*A)
"""
expr = sympify(expr)
is_iterable = iterable(expr)
if not isinstance(expr, Basic) or expr.is_Atom:
if is_iterable:
return type(expr)([factor_terms(i, radical=radical) for i in expr])
return expr
if expr.is_Function or is_iterable or not hasattr(expr, 'args_cnc'):
args = expr.args
newargs = tuple([factor_terms(i, radical=radical) for i in args])
if newargs == args:
return expr
return expr.func(*newargs)
cont, p = expr.as_content_primitive(radical=radical)
list_args, nc = zip(*[ai.args_cnc() for ai in Add.make_args(p)])
list_args = list(list_args)
nc = [((Dummy(), Mul._from_args(i)) if i else None) for i in nc]
ncreps = dict([i for i in nc if i is not None])
for i, a in enumerate(list_args):
if nc[i] is not None:
a.append(nc[i][0])
a = Mul._from_args(a) # gcd_terms will fix up ordering
list_args[i] = gcd_terms(a, isprimitive=True)
# cancel terms that may not have cancelled
p = Add._from_args(list_args) # gcd_terms will fix up ordering
p = gcd_terms(p, isprimitive=True).xreplace(ncreps)
return _keep_coeff(cont, p)
示例2: mask
def mask(terms):
"""replace nc portions of each term with a unique Dummy symbols
and return the replacements to restore them"""
args = [(a, []) if a.is_commutative else a.args_cnc() for a in terms]
reps = []
for i, (c, nc) in enumerate(args):
if nc:
nc = Mul._from_args(nc)
d = Dummy()
reps.append((d, nc))
c.append(d)
args[i] = Mul._from_args(c)
else:
args[i] = c
return args, dict(reps)
示例3: _eval_simplify
def _eval_simplify(self, ratio=1.7, measure=None):
from sympy.simplify.simplify import factor_sum, sum_combine
from sympy.core.function import expand
from sympy.core.mul import Mul
# split the function into adds
terms = Add.make_args(expand(self.function))
s_t = [] # Sum Terms
o_t = [] # Other Terms
for term in terms:
if term.has(Sum):
# if there is an embedded sum here
# it is of the form x * (Sum(whatever))
# hence we make a Mul out of it, and simplify all interior sum terms
subterms = Mul.make_args(expand(term))
out_terms = []
for subterm in subterms:
# go through each term
if isinstance(subterm, Sum):
# if it's a sum, simplify it
out_terms.append(subterm._eval_simplify())
else:
# otherwise, add it as is
out_terms.append(subterm)
# turn it back into a Mul
s_t.append(Mul(*out_terms))
else:
o_t.append(term)
# next try to combine any interior sums for further simplification
result = Add(sum_combine(s_t), *o_t)
return factor_sum(result, limits=self.limits)
示例4: _parse_matrix_expression
def _parse_matrix_expression(expr):
if isinstance(expr, MatMul):
args_nonmat = []
args = []
contractions = []
for arg in expr.args:
if isinstance(arg, MatrixExpr):
args.append(arg)
else:
args_nonmat.append(arg)
contractions = [(2*i+1, 2*i+2) for i in range(len(args)-1)]
return Mul.fromiter(args_nonmat)*CodegenArrayContraction(
CodegenArrayTensorProduct(*[_parse_matrix_expression(arg) for arg in args]),
*contractions
)
elif isinstance(expr, MatAdd):
return CodegenArrayElementwiseAdd(
*[_parse_matrix_expression(arg) for arg in expr.args]
)
elif isinstance(expr, Transpose):
return CodegenArrayPermuteDims(
_parse_matrix_expression(expr.args[0]), [1, 0]
)
else:
return expr
示例5: factor_terms
def factor_terms(expr):
"""Remove common factors from terms in all arguments without
changing the underlying structure of the expr. No expansion or
simplification (and no processing of non-commutative) is performed.
**Examples**
>>> from sympy import factor_terms, Symbol
>>> from sympy.abc import x, y
>>> factor_terms(x + x*(2 + 4*y)**3)
x*(8*(2*y + 1)**3 + 1)
>>> A = Symbol('A', commutative=False)
>>> factor_terms(x*A + x*A + x*y*A)
x*(y*A + 2*A)
"""
expr = sympify(expr)
if iterable(expr):
return type(expr)([factor_terms(i) for i in expr])
if not isinstance(expr, Basic) or expr.is_Atom:
return expr
if expr.is_Function:
return expr.func(*[factor_terms(i) for i in expr.args])
cont, p = expr.as_content_primitive()
list_args, nc = zip(*[ai.args_cnc(clist=True) for ai in Add.make_args(p)])
list_args = list(list_args)
nc = [((Dummy(), Mul._from_args(i)) if i else None) for i in nc]
ncreps = dict([i for i in nc if i is not None])
for i, a in enumerate(list_args):
if nc[i] is not None:
a.append(nc[i][0])
a = Mul._from_args(a) # gcd_terms will fix up ordering
list_args[i] = gcd_terms(a, isprimitive=True)
# cancel terms that may not have cancelled
p = Add._from_args(list_args) # gcd_terms will fix up ordering
p = gcd_terms(p, isprimitive=True).subs(ncreps) # exact subs could be used here
return _keep_coeff(cont, p)
示例6: from_MatMul
def from_MatMul(expr):
args_nonmat = []
args = []
contractions = []
for arg in expr.args:
if isinstance(arg, MatrixExpr):
args.append(arg)
else:
args_nonmat.append(arg)
contractions = [(2*i+1, 2*i+2) for i in range(len(args)-1)]
return Mul.fromiter(args_nonmat)*CodegenArrayContraction(
CodegenArrayTensorProduct(*args),
*contractions
)
示例7: _sqrt
def _sqrt(d):
# remove squares from square root since both will be represented
# in the results; a similar thing is happening in roots() but
# must be duplicated here because not all quadratics are binomials
co = []
other = []
for di in Mul.make_args(d):
if di.is_Pow and di.exp.is_Integer and di.exp % 2 == 0:
co.append(Pow(di.base, di.exp//2))
else:
other.append(di)
if co:
d = Mul(*other)
co = Mul(*co)
return co*sqrt(d)
return sqrt(d)
示例8: eval
def eval(cls, arg):
if arg.is_Number:
if arg is S.NaN:
return S.NaN
elif arg is S.Zero:
return S.One
elif arg is S.One:
return S.Exp1
elif arg is S.Infinity:
return S.Infinity
elif arg is S.NegativeInfinity:
return S.Zero
elif arg.func is log:
return arg.args[0]
elif arg.is_Mul:
Ioo = S.ImaginaryUnit*S.Infinity
if arg in [Ioo, -Ioo]:
return S.NaN
coeff = arg.coeff(S.Pi*S.ImaginaryUnit)
if coeff:
if (2*coeff).is_integer:
if coeff.is_even:
return S.One
elif coeff.is_odd:
return S.NegativeOne
elif (coeff + S.Half).is_even:
return -S.ImaginaryUnit
elif (coeff + S.Half).is_odd:
return S.ImaginaryUnit
# Warning: code in risch.py will be very sensitive to changes
# in this (see DifferentialExtension).
# look for a single log factor
coeff, terms = arg.as_coeff_Mul()
# but it can't be multiplied by oo
if coeff in [S.NegativeInfinity, S.Infinity]:
return None
coeffs, log_term = [coeff], None
for term in Mul.make_args(terms):
if term.func is log:
if log_term is None:
log_term = term.args[0]
else:
return None
elif term.is_comparable:
coeffs.append(term)
else:
return None
return log_term**Mul(*coeffs) if log_term else None
elif arg.is_Add:
out = []
add = []
for a in arg.args:
if a is S.One:
add.append(a)
continue
newa = cls(a)
if newa.func is cls:
add.append(a)
else:
out.append(newa)
if out:
return Mul(*out)*cls(Add(*add), evaluate=False)
elif arg.is_Matrix:
from sympy import Matrix
return arg.exp()
示例9: factor_nc
def factor_nc(expr):
"""Return the factored form of ``expr`` while handling non-commutative
expressions.
**examples**
>>> from sympy.core.exprtools import factor_nc
>>> from sympy import Symbol
>>> from sympy.abc import x
>>> A = Symbol('A', commutative=False)
>>> B = Symbol('B', commutative=False)
>>> factor_nc((x**2 + 2*A*x + A**2).expand())
(x + A)**2
>>> factor_nc(((x + A)*(x + B)).expand())
(x + A)*(x + B)
"""
from sympy.simplify.simplify import _mexpand
from sympy.polys import gcd, factor
expr = sympify(expr)
if not isinstance(expr, Expr) or not expr.args:
return expr
if not expr.is_Add:
return expr.func(*[factor_nc(a) for a in expr.args])
expr, rep, nc_symbols = _mask_nc(expr)
if rep:
return factor(expr).subs(rep)
else:
args = [a.args_cnc() for a in Add.make_args(expr)]
c = g = l = r = S.One
hit = False
# find any commutative gcd term
for i, a in enumerate(args):
if i == 0:
c = Mul._from_args(a[0])
elif a[0]:
c = gcd(c, Mul._from_args(a[0]))
else:
c = S.One
if c is not S.One:
hit = True
c, g = c.as_coeff_Mul()
if g is not S.One:
for i, (cc, _) in enumerate(args):
cc = list(Mul.make_args(Mul._from_args(list(cc))/g))
args[i][0] = cc
else:
for i, (cc, _) in enumerate(args):
cc[0] = cc[0]/c
args[i][0] = cc
# find any noncommutative common prefix
for i, a in enumerate(args):
if i == 0:
n = a[1][:]
else:
n = common_prefix(n, a[1])
if not n:
# is there a power that can be extracted?
if not args[0][1]:
break
b, e = args[0][1][0].as_base_exp()
ok = False
if e.is_Integer:
for t in args:
if not t[1]:
break
bt, et = t[1][0].as_base_exp()
if et.is_Integer and bt == b:
e = min(e, et)
else:
break
else:
ok = hit = True
l = b**e
il = b**-e
for i, a in enumerate(args):
args[i][1][0] = il*args[i][1][0]
break
if not ok:
break
else:
hit = True
lenn = len(n)
l = Mul(*n)
for i, a in enumerate(args):
args[i][1] = args[i][1][lenn:]
# find any noncommutative common suffix
for i, a in enumerate(args):
if i == 0:
n = a[1][:]
else:
n = common_suffix(n, a[1])
if not n:
# is there a power that can be extracted?
if not args[0][1]:
break
b, e = args[0][1][-1].as_base_exp()
ok = False
if e.is_Integer:
for t in args:
#.........这里部分代码省略.........
示例10: __init__
def __init__(self, factors=None): # Factors
"""Initialize Factors from dict or expr.
Examples
========
>>> from sympy.core.exprtools import Factors
>>> from sympy.abc import x
>>> from sympy import I
>>> e = 2*x**3
>>> Factors(e)
Factors({2: 1, x: 3})
>>> Factors(e.as_powers_dict())
Factors({2: 1, x: 3})
>>> f = _
>>> f.factors # underlying dictionary
{2: 1, x: 3}
>>> f.gens # base of each factor
frozenset([2, x])
>>> Factors(0)
Factors({0: 1})
>>> Factors(I)
Factors({I: 1})
Notes
=====
Although a dictionary can be passed, only minimal checking is
performed: powers of -1 and I are made canonical.
"""
if isinstance(factors, (SYMPY_INTS, float)):
factors = S(factors)
if isinstance(factors, Factors):
factors = factors.factors.copy()
elif factors is None or factors is S.One:
factors = {}
elif factors is S.Zero or factors == 0:
factors = {S.Zero: S.One}
elif isinstance(factors, Number):
n = factors
factors = {}
if n < 0:
factors[S.NegativeOne] = S.One
n = -n
if n is not S.One:
if n.is_Float or n.is_Integer or n is S.Infinity:
factors[n] = S.One
elif n.is_Rational:
# since we're processing Numbers, the denominator is
# stored with a negative exponent; all other factors
# are left .
if n.p != 1:
factors[Integer(n.p)] = S.One
factors[Integer(n.q)] = S.NegativeOne
else:
raise ValueError('Expected Float|Rational|Integer, not %s' % n)
elif isinstance(factors, Basic) and not factors.args:
factors = {factors: S.One}
elif isinstance(factors, Expr):
c, nc = factors.args_cnc()
i = c.count(I)
for _ in range(i):
c.remove(I)
factors = dict(Mul._from_args(c).as_powers_dict())
if i:
factors[I] = S.One*i
if nc:
factors[Mul(*nc, evaluate=False)] = S.One
else:
factors = factors.copy() # /!\ should be dict-like
# tidy up -/+1 and I exponents if Rational
handle = []
for k in factors:
if k is I or k in (-1, 1):
handle.append(k)
if handle:
i1 = S.One
for k in handle:
if not _isnumber(factors[k]):
continue
i1 *= k**factors.pop(k)
if i1 is not S.One:
for a in i1.args if i1.is_Mul else [i1]: # at worst, -1.0*I*(-1)**e
if a is S.NegativeOne:
factors[a] = S.One
elif a is I:
factors[I] = S.One
elif a.is_Pow:
if S.NegativeOne not in factors:
factors[S.NegativeOne] = S.Zero
factors[S.NegativeOne] += a.exp
elif a == 1:
factors[a] = S.One
elif a == -1:
factors[-a] = S.One
factors[S.NegativeOne] = S.One
#.........这里部分代码省略.........
示例11: rsolve
def rsolve(f, y, init=None):
"""
Solve univariate recurrence with rational coefficients.
Given `k`-th order linear recurrence `\operatorname{L} y = f`,
or equivalently:
.. math:: a_{k}(n) y(n+k) + a_{k-1}(n) y(n+k-1) +
\dots + a_{0}(n) y(n) = f(n)
where `a_{i}(n)`, for `i=0, \dots, k`, are polynomials or rational
functions in `n`, and `f` is a hypergeometric function or a sum
of a fixed number of pairwise dissimilar hypergeometric terms in
`n`, finds all solutions or returns ``None``, if none were found.
Initial conditions can be given as a dictionary in two forms:
(1) ``{ n_0 : v_0, n_1 : v_1, ..., n_m : v_m }``
(2) ``{ y(n_0) : v_0, y(n_1) : v_1, ..., y(n_m) : v_m }``
or as a list ``L`` of values:
``L = [ v_0, v_1, ..., v_m ]``
where ``L[i] = v_i``, for `i=0, \dots, m`, maps to `y(n_i)`.
Examples
========
Lets consider the following recurrence:
.. math:: (n - 1) y(n + 2) - (n^2 + 3 n - 2) y(n + 1) +
2 n (n + 1) y(n) = 0
>>> from sympy import Function, rsolve
>>> from sympy.abc import n
>>> y = Function('y')
>>> f = (n - 1)*y(n + 2) - (n**2 + 3*n - 2)*y(n + 1) + 2*n*(n + 1)*y(n)
>>> rsolve(f, y(n))
2**n*C0 + C1*factorial(n)
>>> rsolve(f, y(n), { y(0):0, y(1):3 })
3*2**n - 3*factorial(n)
See Also
========
rsolve_poly, rsolve_ratio, rsolve_hyper
"""
if isinstance(f, Equality):
f = f.lhs - f.rhs
n = y.args[0]
k = Wild('k', exclude=(n,))
# Preprocess user input to allow things like
# y(n) + a*(y(n + 1) + y(n - 1))/2
f = f.expand().collect(y.func(Wild('m', integer=True)))
h_part = defaultdict(lambda: S.Zero)
i_part = S.Zero
for g in Add.make_args(f):
coeff = S.One
kspec = None
for h in Mul.make_args(g):
if h.is_Function:
if h.func == y.func:
result = h.args[0].match(n + k)
if result is not None:
kspec = int(result[k])
else:
raise ValueError(
"'%s(%s+k)' expected, got '%s'" % (y.func, n, h))
else:
raise ValueError(
"'%s' expected, got '%s'" % (y.func, h.func))
else:
coeff *= h
if kspec is not None:
h_part[kspec] += coeff
else:
i_part += coeff
for k, coeff in h_part.iteritems():
h_part[k] = simplify(coeff)
common = S.One
for coeff in h_part.itervalues():
if coeff.is_rational_function(n):
if not coeff.is_polynomial(n):
common = lcm(common, coeff.as_numer_denom()[1], n)
else:
raise ValueError(
"Polynomial or rational function expected, got '%s'" % coeff)
#.........这里部分代码省略.........
示例12: eval
def eval(cls, p, q):
from sympy.core.add import Add
from sympy.core.mul import Mul
from sympy.core.singleton import S
from sympy.core.exprtools import gcd_terms
from sympy.polys.polytools import gcd
def doit(p, q):
"""Try to return p % q if both are numbers or +/-p is known
to be less than or equal q.
"""
if q == S.Zero:
raise ZeroDivisionError("Modulo by zero")
if p.is_infinite or q.is_infinite or p is nan or q is nan:
return nan
if p == S.Zero or p == q or p == -q or (p.is_integer and q == 1):
return S.Zero
if q.is_Number:
if p.is_Number:
return p%q
if q == 2:
if p.is_even:
return S.Zero
elif p.is_odd:
return S.One
if hasattr(p, '_eval_Mod'):
rv = getattr(p, '_eval_Mod')(q)
if rv is not None:
return rv
# by ratio
r = p/q
try:
d = int(r)
except TypeError:
pass
else:
if isinstance(d, integer_types):
rv = p - d*q
if (rv*q < 0) == True:
rv += q
return rv
# by difference
# -2|q| < p < 2|q|
d = abs(p)
for _ in range(2):
d -= abs(q)
if d.is_negative:
if q.is_positive:
if p.is_positive:
return d + q
elif p.is_negative:
return -d
elif q.is_negative:
if p.is_positive:
return d
elif p.is_negative:
return -d + q
break
rv = doit(p, q)
if rv is not None:
return rv
# denest
if isinstance(p, cls):
qinner = p.args[1]
if qinner % q == 0:
return cls(p.args[0], q)
elif (qinner*(q - qinner)).is_nonnegative:
# |qinner| < |q| and have same sign
return p
elif isinstance(-p, cls):
qinner = (-p).args[1]
if qinner % q == 0:
return cls(-(-p).args[0], q)
elif (qinner*(q + qinner)).is_nonpositive:
# |qinner| < |q| and have different sign
return p
elif isinstance(p, Add):
# separating into modulus and non modulus
both_l = non_mod_l, mod_l = [], []
for arg in p.args:
both_l[isinstance(arg, cls)].append(arg)
# if q same for all
if mod_l and all(inner.args[1] == q for inner in mod_l):
net = Add(*non_mod_l) + Add(*[i.args[0] for i in mod_l])
return cls(net, q)
elif isinstance(p, Mul):
# separating into modulus and non modulus
both_l = non_mod_l, mod_l = [], []
for arg in p.args:
both_l[isinstance(arg, cls)].append(arg)
if mod_l and all(inner.args[1] == q for inner in mod_l):
#.........这里部分代码省略.........
示例13: _monotonic_sign
#.........这里部分代码省略.........
from sympy.polys.polytools import real_roots
from sympy.polys.polyroots import roots
from sympy.polys.polyerrors import PolynomialError
x = free.pop()
x0 = _monotonic_sign(x)
if x0 == _eps or x0 == -_eps:
x0 = S.Zero
if x0 is not None:
d = self.diff(x)
if d.is_number:
roots = []
else:
try:
roots = real_roots(d)
except (PolynomialError, NotImplementedError):
roots = [r for r in roots(d, x) if r.is_real]
y = self.subs(x, x0)
if x.is_nonnegative and all(r <= x0 for r in roots):
if y.is_nonnegative and d.is_positive:
if y:
return y if y.is_positive else Dummy('pos', positive=True)
else:
return Dummy('nneg', nonnegative=True)
if y.is_nonpositive and d.is_negative:
if y:
return y if y.is_negative else Dummy('neg', negative=True)
else:
return Dummy('npos', nonpositive=True)
elif x.is_nonpositive and all(r >= x0 for r in roots):
if y.is_nonnegative and d.is_negative:
if y:
return Dummy('pos', positive=True)
else:
return Dummy('nneg', nonnegative=True)
if y.is_nonpositive and d.is_positive:
if y:
return Dummy('neg', negative=True)
else:
return Dummy('npos', nonpositive=True)
else:
n, d = self.as_numer_denom()
den = None
if n.is_number:
den = _monotonic_sign(d)
elif not d.is_number:
if _monotonic_sign(n) is not None:
den = _monotonic_sign(d)
if den is not None and (den.is_positive or den.is_negative):
v = n*den
if v.is_positive:
return Dummy('pos', positive=True)
elif v.is_nonnegative:
return Dummy('nneg', nonnegative=True)
elif v.is_negative:
return Dummy('neg', negative=True)
elif v.is_nonpositive:
return Dummy('npos', nonpositive=True)
return None
# multivariate
c, a = self.as_coeff_Add()
v = None
if not a.is_polynomial():
# F/A or A/F where A is a number and F is a signed, rational monomial
n, d = a.as_numer_denom()
if not (n.is_number or d.is_number):
return
if (
a.is_Mul or a.is_Pow) and \
a.is_rational and \
all(p.exp.is_Integer for p in a.atoms(Pow) if p.is_Pow) and \
(a.is_positive or a.is_negative):
v = S(1)
for ai in Mul.make_args(a):
if ai.is_number:
v *= ai
continue
reps = {}
for x in ai.free_symbols:
reps[x] = _monotonic_sign(x)
if reps[x] is None:
return
v *= ai.subs(reps)
elif c:
# signed linear expression
if not any(p for p in a.atoms(Pow) if not p.is_number) and (a.is_nonpositive or a.is_nonnegative):
free = list(a.free_symbols)
p = {}
for i in free:
v = _monotonic_sign(i)
if v is None:
return
p[i] = v or (_eps if i.is_nonnegative else -_eps)
v = a.xreplace(p)
if v is not None:
rv = v + c
if v.is_nonnegative and rv.is_positive:
return rv.subs(_eps, 0)
if v.is_nonpositive and rv.is_negative:
return rv.subs(_eps, 0)
示例14: eval
def eval(cls, arg):
from sympy.assumptions import ask, Q
from sympy.calculus import AccumBounds
if arg.is_Number:
if arg is S.NaN:
return S.NaN
elif arg is S.Zero:
return S.One
elif arg is S.One:
return S.Exp1
elif arg is S.Infinity:
return S.Infinity
elif arg is S.NegativeInfinity:
return S.Zero
elif isinstance(arg, log):
return arg.args[0]
elif isinstance(arg, AccumBounds):
return AccumBounds(exp(arg.min), exp(arg.max))
elif arg.is_Mul:
if arg.is_number or arg.is_Symbol:
coeff = arg.coeff(S.Pi*S.ImaginaryUnit)
if coeff:
if ask(Q.integer(2*coeff)):
if ask(Q.even(coeff)):
return S.One
elif ask(Q.odd(coeff)):
return S.NegativeOne
elif ask(Q.even(coeff + S.Half)):
return -S.ImaginaryUnit
elif ask(Q.odd(coeff + S.Half)):
return S.ImaginaryUnit
# Warning: code in risch.py will be very sensitive to changes
# in this (see DifferentialExtension).
# look for a single log factor
coeff, terms = arg.as_coeff_Mul()
# but it can't be multiplied by oo
if coeff in [S.NegativeInfinity, S.Infinity]:
return None
coeffs, log_term = [coeff], None
for term in Mul.make_args(terms):
if isinstance(term, log):
if log_term is None:
log_term = term.args[0]
else:
return None
elif term.is_comparable:
coeffs.append(term)
else:
return None
return log_term**Mul(*coeffs) if log_term else None
elif arg.is_Add:
out = []
add = []
for a in arg.args:
if a is S.One:
add.append(a)
continue
newa = cls(a)
if isinstance(newa, cls):
add.append(a)
else:
out.append(newa)
if out:
return Mul(*out)*cls(Add(*add), evaluate=False)
elif arg.is_Matrix:
return arg.exp()
示例15: rsolve
def rsolve(f, y, init=None):
"""
Solve univariate recurrence with rational coefficients.
Given k-th order linear recurrence Ly = f, or equivalently:
a_{k}(n) y(n+k) + a_{k-1}(n) y(n+k-1) + ... + a_{0}(n) y(n) = f
where a_{i}(n), for i=0..k, are polynomials or rational functions
in n, and f is a hypergeometric function or a sum of a fixed number
of pairwise dissimilar hypergeometric terms in n, finds all solutions
or returns None, if none were found.
Initial conditions can be given as a dictionary in two forms:
[1] { n_0 : v_0, n_1 : v_1, ..., n_m : v_m }
[2] { y(n_0) : v_0, y(n_1) : v_1, ..., y(n_m) : v_m }
or as a list L of values:
L = [ v_0, v_1, ..., v_m ]
where L[i] = v_i, for i=0..m, maps to y(n_i).
As an example lets consider the following recurrence:
(n - 1) y(n + 2) - (n**2 + 3 n - 2) y(n + 1) + 2 n (n + 1) y(n) == 0
>>> from sympy import Function, rsolve
>>> from sympy.abc import n
>>> y = Function('y')
>>> f = (n-1)*y(n+2) - (n**2+3*n-2)*y(n+1) + 2*n*(n+1)*y(n)
>>> rsolve(f, y(n))
2**n*C0 + C1*n!
>>> rsolve(f, y(n), { y(0):0, y(1):3 })
3*2**n - 3*n!
"""
if isinstance(f, Equality):
f = f.lhs - f.rhs
n = y.args[0]
k = Wild('k', exclude=(n,))
h_part = defaultdict(lambda: S.Zero)
i_part = S.Zero
for g in Add.make_args(f):
coeff = S.One
kspec = None
for h in Mul.make_args(g):
if h.is_Function:
if h.func == y.func:
result = h.args[0].match(n + k)
if result is not None:
kspec = int(result[k])
else:
raise ValueError(
"'%s(%s+k)' expected, got '%s'" % (y.func, n, h))
else:
raise ValueError(
"'%s' expected, got '%s'" % (y.func, h.func))
else:
coeff *= h
if kspec is not None:
h_part[kspec] += coeff
else:
i_part += coeff
for k, coeff in h_part.iteritems():
h_part[k] = simplify(coeff)
common = S.One
for coeff in h_part.itervalues():
if coeff.is_rational_function(n):
if not coeff.is_polynomial(n):
common = lcm(common, coeff.as_numer_denom()[1], n)
else:
raise ValueError(
"Polynomial or rational function expected, got '%s'" % coeff)
i_numer, i_denom = i_part.as_numer_denom()
if i_denom.is_polynomial(n):
common = lcm(common, i_denom, n)
if common is not S.One:
for k, coeff in h_part.iteritems():
numer, denom = coeff.as_numer_denom()
h_part[k] = numer*quo(common, denom, n)
i_part = i_numer*quo(common, i_denom, n)
K_min = min(h_part.keys())
#.........这里部分代码省略.........