sympy简明用法

摘要:
系统地学习什么是SympySympy是一个可以执行符号运算的第三方科学计算库。数学对象可以精确地表示而不是近似值,这也意味着没有计算的未知量可以以符号的形式留在数学表达式中。Sympy.sqrt#更好的打印效果[sqrt{3}]一个更有趣的例子是使用符号来定义变量,也就是说,您必须在使用变量之前定义它。
系统学习Sympy

什么是Sympy

Sympy 是一个可以进行符号运算的第三方科学计算库,数学对象可以被精确的表达,而不是近似值,这也意味着带有未计算的未知量可以以符号的形式留在数学表达式中。

import sympy
sympy.sqrt(3) #用Sympy表达无理数
sqrt(3)
sympy.init_printing(use_unicode=True) # `use_unicode=True` 是默认值,在程序最前面调用该函数,可以用Unicode将结果打印的更好看。
sympy.sqrt(3) #更好的print 效果

[sqrt{3} ]

一个更有趣的例子

在Sympy 中,用symbols来定义变量,也就是说,在使用某个变量前,必须先定义它。

x,y=sympy.symbols('x,y') # symbols('x y') 或者symbols('x,y')
help(sympy.symbols)
Help on function symbols in module sympy.core.symbol:

symbols(names, **args)
    Transform strings into instances of :class:`Symbol` class.
    
    :func:`symbols` function returns a sequence of symbols with names taken
    from ``names`` argument, which can be a comma or whitespace delimited
    string, or a sequence of strings::
    
        >>> from sympy import symbols, Function
    
        >>> x, y, z = symbols('x,y,z')
        >>> a, b, c = symbols('a b c')
    
    The type of output is dependent on the properties of input arguments::
    
        >>> symbols('x')
        x
        >>> symbols('x,')
        (x,)
        >>> symbols('x,y')
        (x, y)
        >>> symbols(('a', 'b', 'c'))
        (a, b, c)
        >>> symbols(['a', 'b', 'c'])
        [a, b, c]
        >>> symbols(set(['a', 'b', 'c']))
        set([a, b, c])
    
    If an iterable container is needed for a single symbol, set the ``seq``
    argument to ``True`` or terminate the symbol name with a comma::
    
        >>> symbols('x', seq=True)
        (x,)
    
    To reduce typing, range syntax is supported to create indexed symbols.
    Ranges are indicated by a colon and the type of range is determined by
    the character to the right of the colon. If the character is a digit
    then all contiguous digits to the left are taken as the nonnegative
    starting value (or 0 if there is no digit left of the colon) and all
    contiguous digits to the right are taken as 1 greater than the ending
    value::
    
        >>> symbols('x:10')
        (x0, x1, x2, x3, x4, x5, x6, x7, x8, x9)
    
        >>> symbols('x5:10')
        (x5, x6, x7, x8, x9)
        >>> symbols('x5(:2)')
        (x50, x51)
    
        >>> symbols('x5:10,y:5')
        (x5, x6, x7, x8, x9, y0, y1, y2, y3, y4)
    
        >>> symbols(('x5:10', 'y:5'))
        ((x5, x6, x7, x8, x9), (y0, y1, y2, y3, y4))
    
    If the character to the right of the colon is a letter, then the single
    letter to the left (or 'a' if there is none) is taken as the start
    and all characters in the lexicographic range *through* the letter to
    the right are used as the range::
    
        >>> symbols('x:z')
        (x, y, z)
        >>> symbols('x:c')  # null range
        ()
        >>> symbols('x(:c)')
        (xa, xb, xc)
    
        >>> symbols(':c')
        (a, b, c)
    
        >>> symbols('a:d, x:z')
        (a, b, c, d, x, y, z)
    
        >>> symbols(('a:d', 'x:z'))
        ((a, b, c, d), (x, y, z))
    
    Multiple ranges are supported; contiguous numerical ranges should be
    separated by parentheses to disambiguate the ending number of one
    range from the starting number of the next::
    
        >>> symbols('x:2(1:3)')
        (x01, x02, x11, x12)
        >>> symbols(':3:2')  # parsing is from left to right
        (00, 01, 10, 11, 20, 21)
    
    Only one pair of parentheses surrounding ranges are removed, so to
    include parentheses around ranges, double them. And to include spaces,
    commas, or colons, escape them with a backslash::
    
        >>> symbols('x((a:b))')
        (x(a), x(b))
        >>> symbols('x(:1\,:2)')  # or 'x((:1)\,(:2))'
        (x(0,0), x(0,1))
    
    All newly created symbols have assumptions set according to ``args``::
    
        >>> a = symbols('a', integer=True)
        >>> a.is_integer
        True
    
        >>> x, y, z = symbols('x,y,z', real=True)
        >>> x.is_real and y.is_real and z.is_real
        True
    
    Despite its name, :func:`symbols` can create symbol-like objects like
    instances of Function or Wild classes. To achieve this, set ``cls``
    keyword argument to the desired type::
    
        >>> symbols('f,g,h', cls=Function)
        (f, g, h)
    
        >>> type(_[0])
        <class 'sympy.core.function.UndefinedFunction'>
x,y=sympy.symbols('x y')
expr=x+2*y;expr

[x + 2 y ]

expr+1

[x + 2 y + 1 ]

x*expr # 除了 x-x=0,sqr(8)=2*sqrt(2),很多简化不能自动的进行,这是因为我们有可能希望展开,也有可能希望是多因式的形式,在sympy中,有专门的函数来进行转化

[x left(x + 2 y ight) ]

expanded_expr=sympy.expand(x*expr)
expanded_expr

[x^{2} + 2 x y ]

sympy.factor(expanded_expr)

[x left(x + 2 y ight) ]

强大的Sympy

Sympy的强大在于可以进行很多种符号运算,比如简化表达式,计算导数,积分,极限,解方程,矩阵运算。下面简单的给几个例子:

x,y,z,nu=sympy.symbols('x y z nu')
  • 对 sin(x)exp(x)求导
sympy.diff(sympy.sin(x)*sympy.exp(x),x)

[e^{x} sin{left (x ight )} + e^{x} cos{left (x ight )} ]

  • 计算exsin(x)+excos(x)的积分
sympy.integrate(sympy.exp(x)*(sympy.sin(x)+sympy.cos(x)),x)

[e^{x} sin{left (x ight )} ]

  • 计算sin(x2)的无穷积分:
from sympy import *
sympy.integrate(sympy.sin(x**2),(x,-oo,oo))

[frac{sqrt{2} sqrt{pi}}{2} ]

  • sin(x)/x,当x趋近于0时的极限:
sympy.limit(sympy.sin(x)/x,x,0)

[1 ]

  • 解方程:x**2-2=0
sympy.solve(x**2-2,x)

[left [ - sqrt{2}, quad sqrt{2} ight ] ]

  • 解微分方程 y''-y=e(t)
y=sympy.Function('y')
t=sympy.symbols('t')
sympy.dsolve(sympy.Eq(y(t).diff(t,t)-y(t),sympy.exp(t)),y(t))

[y{left (t ight )} = C_{2} e^{- t} + left(C_{1} + frac{t}{2} ight) e^{t} ]

  • 找矩阵[[1,2],[2,2]]的特征值:
sympy.Matrix([[1,2],[2,2]]).eigenvals()

[left { frac{3}{2} + frac{sqrt{17}}{2} : 1, quad - frac{sqrt{17}}{2} + frac{3}{2} : 1 ight } ]

  • 用LATEX打印cos(x)**2在(0,pi)的积分:
sympy.integrate(sympy.cos(x)**2,(x,0,sympy.pi))

[frac{pi}{2} ]

sympy.integrate(sympy.cos(x)**2,x)

[frac{x}{2} + frac{1}{2} sin{left (x ight )} cos{left (x ight )} ]

sympy.latex(sympy.Integral(sympy.cos(x)**2,(x,0,sympy.pi)))
'\int_{0}^{\pi} \cos^{2}{\left (x \right )}\, dx'
sympy.Integral(sympy.cos(x)**2,(x,0,sympy.pi))

[int_{0}^{pi} cos^{2}{left (x ight )}\, dx ]

Attention:Integral is used for the showing purpose,while integrate is used for evaluating purpose.

下面通过一个小例子来说明在Sympy中的symbols和Python中的变量的区别:

x=sympy.symbols('x')
expr=x+1
x=2
# expr
type(expr)
sympy.core.add.Add

可以发现expr并不是3,因为x=2只是改变了x的指向,并不影响到expr,如果要想将expr中的x替换为2,可以用expr的subs方法。

expr.subs(x,2)

[x + 1 ]

等号

在Python中,=代表赋值,为了不与python冲突,所以等式不能用=来表示,而==代表逻辑判断等于,其返回布尔值,显然作为表示符号运算中的“等于”,也不合适。在Sympy中,我们用Eq来表示等式:

sympy.Eq(x+1,4)

[mathrm{False} ]

为了比较两个式子是否相等,可以将两个式子做减法,然后用simplify方法:

x=sympy.symbols('x')
a=(x+1)**2;a

[left(x + 1 ight)^{2} ]

b=x**2+1+2*x;b

[x^{2} + 2 x + 1 ]

(a-b).simplify()

[0 ]

另外,还有一个equals方法可以用来判断,两个式子是否是相同的:

a.equals(b)
True

^和/

我们用**来表示指数,而不是^,这是因为Sympy跟Python保持一致,在Python中,^表示逻辑运算符‘异或’:

True ^ False
True
True ^ True
False
False ^ False
False

基本操作

这里先简单介绍下Sympy中的基本操作:

x,y,z=sympy.symbols('x y z')

替换

help(expr.subs)
Help on method subs in module sympy.core.basic:

subs(*args, **kwargs) method of sympy.core.add.Add instance
    Substitutes old for new in an expression after sympifying args.
    
    `args` is either:
      - two arguments, e.g. foo.subs(old, new)
      - one iterable argument, e.g. foo.subs(iterable). The iterable may be
         o an iterable container with (old, new) pairs. In this case the
           replacements are processed in the order given with successive
           patterns possibly affecting replacements already made.
         o a dict or set whose key/value items correspond to old/new pairs.
           In this case the old/new pairs will be sorted by op count and in
           case of a tie, by number of args and the default_sort_key. The
           resulting sorted list is then processed as an iterable container
           (see previous).
    
    If the keyword ``simultaneous`` is True, the subexpressions will not be
    evaluated until all the substitutions have been made.
    
    Examples
    ========
    
    >>> from sympy import pi, exp, limit, oo
    >>> from sympy.abc import x, y
    >>> (1 + x*y).subs(x, pi)
    pi*y + 1
    >>> (1 + x*y).subs({x:pi, y:2})
    1 + 2*pi
    >>> (1 + x*y).subs([(x, pi), (y, 2)])
    1 + 2*pi
    >>> reps = [(y, x**2), (x, 2)]
    >>> (x + y).subs(reps)
    6
    >>> (x + y).subs(reversed(reps))
    x**2 + 2
    
    >>> (x**2 + x**4).subs(x**2, y)
    y**2 + y
    
    To replace only the x**2 but not the x**4, use xreplace:
    
    >>> (x**2 + x**4).xreplace({x**2: y})
    x**4 + y
    
    To delay evaluation until all substitutions have been made,
    set the keyword ``simultaneous`` to True:
    
    >>> (x/y).subs([(x, 0), (y, 0)])
    0
    >>> (x/y).subs([(x, 0), (y, 0)], simultaneous=True)
    nan
    
    This has the added feature of not allowing subsequent substitutions
    to affect those already made:
    
    >>> ((x + y)/y).subs({x + y: y, y: x + y})
    1
    >>> ((x + y)/y).subs({x + y: y, y: x + y}, simultaneous=True)
    y/(x + y)
    
    In order to obtain a canonical result, unordered iterables are
    sorted by count_op length, number of arguments and by the
    default_sort_key to break any ties. All other iterables are left
    unsorted.
    
    >>> from sympy import sqrt, sin, cos
    >>> from sympy.abc import a, b, c, d, e
    
    >>> A = (sqrt(sin(2*x)), a)
    >>> B = (sin(2*x), b)
    >>> C = (cos(2*x), c)
    >>> D = (x, d)
    >>> E = (exp(x), e)
    
    >>> expr = sqrt(sin(2*x))*sin(exp(x)*x)*cos(2*x) + sin(2*x)
    
    >>> expr.subs(dict([A, B, C, D, E]))
    a*c*sin(d*e) + b
    
    The resulting expression represents a literal replacement of the
    old arguments with the new arguments. This may not reflect the
    limiting behavior of the expression:
    
    >>> (x**3 - 3*x).subs({x: oo})
    nan
    
    >>> limit(x**3 - 3*x, x, oo)
    oo
    
    If the substitution will be followed by numerical
    evaluation, it is better to pass the substitution to
    evalf as
    
    >>> (1/x).evalf(subs={x: 3.0}, n=21)
    0.333333333333333333333
    
    rather than
    
    >>> (1/x).subs({x: 3.0}).evalf(21)
    0.333333333333333314830
    
    as the former will ensure that the desired level of precision is
    obtained.
    
    See Also
    ========
    replace: replacement capable of doing wildcard-like matching,
             parsing of match, and conditional replacements
    xreplace: exact node replacement in expr tree; also capable of
              using matching rules
    evalf: calculates the given formula to a desired level of precision
expr=sympy.cos(x)+1;expr

[cos{left (x ight )} + 1 ]

expr.subs(x,y)

[cos{left (y ight )} + 1 ]

subs主要用来做以下两个事:

  • 计算表达式在某个值的结果.
(x+y).subs([(x,1),(y,2)])

[3 ]

  • 用子表达式来代替另一个子表达式
expr=x**y;expr

[x^{y} ]

expr=expr.subs(y,x**x);expr

[x^{x^{x}} ]

expr=sympy.sin(2*x)+sympy.cos(2*x);expr

[sin{left (2 x ight )} + cos{left (2 x ight )} ]

sympy.expand_trig(expr)

[2 sin{left (x ight )} cos{left (x ight )} + 2 cos^{2}{left (x ight )} - 1 ]

expr.subs(sympy.sin(2*x),sympy.sin(x)*sympy.cos(x)*2)

[2 sin{left (x ight )} cos{left (x ight )} + cos{left (2 x ight )} ]

有两点注意事项:(1)subs返回一个新的表达式,sympy对象都是不可变的,这意味着subs没有在原处修改:

expr=sympy.cos(x)
expr.subs(x,0)

[1 ]

(2)传入list(old1,new1),(old2,new2)来一次性完成多个替换:

(x**3+4*x*y-z).subs([(x,2),(y,4),(z,0)])

[40 ]

将字符串转化为Sympy表达式

simpify函数(注意不是simplify) 用来将字符串转化为Sympy中的表达式:

str_expr='x**2+3*x-1/2'
sympy.sympify(str_expr)

[x^{2} + 3 x - frac{1}{2} ]

sympy.sympify(str_expr).subs(x,2)

[frac{19}{2} ]

注意 sympify 用了eval,对于没初始化的变量,不要用

sympy.sympify('xxx+xyy')

[xxx + xyy ]

# sympy.sympify('xxx+xyy').subs(xxx,3) 

evalf

evalf主要用来把数值表达式来计算为浮点数

sympy.sqrt(8).evalf()

[2.82842712474619 ]

sympy.sqrt(8).evalf(3) #结果保留3位有效数字

[2.83 ]

expr=sympy.cos(2*x);expr

[cos{left (2 x ight )} ]

expr.subs(x,2)

[cos{left (4 ight )} ]

expr.subs(x,2).evalf()

[-0.653643620863612 ]

求解一个数值表达式的某点的具体值,一般是先subs,然后evalf,但还有更高效的方法:在evalf中,传入subs=dict(),其中dict是一对symbol:point.

expr

[cos{left (2 x ight )} ]

expr.evalf(subs={x:2})

[-0.653643620863612 ]

有时候,由于在计算完成后,存在小于一定程度的精度的四舍五入误差,可以通过设置chop为True:

one=sympy.sin(1)**2+sympy.cos(1)**2
one

[cos^{2}{left (1 ight )} + sin^{2}{left (1 ight )} ]

(one-1).evalf(chop=True)

[0 ]

(one-1).evalf()

[-4.0 cdot 10^{-124} ]

lambdify

subsevalf对于计算表达式在一个点的值比较方便,但如果想同时求很多个点,用sympy就会非常慢,这种情况下,使用其他库,如numpy,scipy等就比较方便,而将sympy表达式转换为函数的就是lambdify方法,lambdify就像lambda,除了它还需要将sympy中的变量转换为给定库中的变量:

import numpy as np
a=np.arange(10)
expr=sympy.sin(x)
f=sympy.lambdify(x,expr,'numpy')
f(a)
array([ 0.        ,  0.84147098,  0.90929743,  0.14112001, -0.7568025 ,
       -0.95892427, -0.2794155 ,  0.6569866 ,  0.98935825,  0.41211849])
help(sympy.lambdify)
Help on function lambdify in module sympy.utilities.lambdify:

lambdify(args, expr, modules=None, printer=None, use_imps=True, dummify=True)
    Returns a lambda function for fast calculation of numerical values.
    
    If not specified differently by the user, SymPy functions are replaced as
    far as possible by either python-math, numpy (if available) or mpmath
    functions - exactly in this order. To change this behavior, the "modules"
    argument can be used. It accepts:
    
     - the strings "math", "mpmath", "numpy", "numexpr", "sympy"
     - any modules (e.g. math)
     - dictionaries that map names of sympy functions to arbitrary functions
     - lists that contain a mix of the arguments above, with higher priority
       given to entries appearing first.
    
    The default behavior is to substitute all arguments in the provided
    expression with dummy symbols. This allows for applied functions (e.g.
    f(t)) to be supplied as arguments. Call the function with dummify=False if
    dummy substitution is unwanted (and `args` is not a string). If you want
    to view the lambdified function or provide "sympy" as the module, you
    should probably set dummify=False.
    
    For functions involving large array calculations, numexpr can provide a
    significant speedup over numpy.  Please note that the available functions
    for numexpr are more limited than numpy but can be expanded with
    implemented_function and user defined subclasses of Function.  If specified,
    numexpr may be the only option in modules. The official list of numexpr
    functions can be found at:
    https://github.com/pydata/numexpr#supported-functions
    
    In previous releases ``lambdify`` replaced ``Matrix`` with ``numpy.matrix``
    by default. As of release 1.0 ``numpy.array`` is the default.
    To get the old default behavior you must pass in ``[{'ImmutableMatrix':
    numpy.matrix}, 'numpy']`` to the ``modules`` kwarg.
    
    >>> from sympy import lambdify, Matrix
    >>> from sympy.abc import x, y
    >>> import numpy
    >>> array2mat = [{'ImmutableMatrix': numpy.matrix}, 'numpy']
    >>> f = lambdify((x, y), Matrix([x, y]), modules=array2mat)
    >>> f(1, 2)
    matrix([[1],
            [2]])
    
    Usage
    =====
    
    (1) Use one of the provided modules:
    
        >>> from sympy import sin, tan, gamma
        >>> from sympy.utilities.lambdify import lambdastr
        >>> from sympy.abc import x, y
        >>> f = lambdify(x, sin(x), "math")
    
        Attention: Functions that are not in the math module will throw a name
                   error when the lambda function is evaluated! So this would
                   be better:
    
        >>> f = lambdify(x, sin(x)*gamma(x), ("math", "mpmath", "sympy"))
    
    (2) Use some other module:
    
        >>> import numpy
        >>> f = lambdify((x,y), tan(x*y), numpy)
    
        Attention: There are naming differences between numpy and sympy. So if
                   you simply take the numpy module, e.g. sympy.atan will not be
                   translated to numpy.arctan. Use the modified module instead
                   by passing the string "numpy":
    
        >>> f = lambdify((x,y), tan(x*y), "numpy")
        >>> f(1, 2)
        -2.18503986326
        >>> from numpy import array
        >>> f(array([1, 2, 3]), array([2, 3, 5]))
        [-2.18503986 -0.29100619 -0.8559934 ]
    
    (3) Use a dictionary defining custom functions:
    
        >>> def my_cool_function(x): return 'sin(%s) is cool' % x
        >>> myfuncs = {"sin" : my_cool_function}
        >>> f = lambdify(x, sin(x), myfuncs); f(1)
        'sin(1) is cool'
    
    Examples
    ========
    
    >>> from sympy.utilities.lambdify import implemented_function
    >>> from sympy import sqrt, sin, Matrix
    >>> from sympy import Function
    >>> from sympy.abc import w, x, y, z
    
    >>> f = lambdify(x, x**2)
    >>> f(2)
    4
    >>> f = lambdify((x, y, z), [z, y, x])
    >>> f(1,2,3)
    [3, 2, 1]
    >>> f = lambdify(x, sqrt(x))
    >>> f(4)
    2.0
    >>> f = lambdify((x, y), sin(x*y)**2)
    >>> f(0, 5)
    0.0
    >>> row = lambdify((x, y), Matrix((x, x + y)).T, modules='sympy')
    >>> row(1, 2)
    Matrix([[1, 3]])
    
    Tuple arguments are handled and the lambdified function should
    be called with the same type of arguments as were used to create
    the function.:
    
    >>> f = lambdify((x, (y, z)), x + y)
    >>> f(1, (2, 4))
    3
    
    A more robust way of handling this is to always work with flattened
    arguments:
    
    >>> from sympy.utilities.iterables import flatten
    >>> args = w, (x, (y, z))
    >>> vals = 1, (2, (3, 4))
    >>> f = lambdify(flatten(args), w + x + y + z)
    >>> f(*flatten(vals))
    10
    
    Functions present in `expr` can also carry their own numerical
    implementations, in a callable attached to the ``_imp_``
    attribute.  Usually you attach this using the
    ``implemented_function`` factory:
    
    >>> f = implemented_function(Function('f'), lambda x: x+1)
    >>> func = lambdify(x, f(x))
    >>> func(4)
    5
    
    ``lambdify`` always prefers ``_imp_`` implementations to implementations
    in other namespaces, unless the ``use_imps`` input parameter is False.

简化

simplify

sympy有很多方法来进行简化操作,也有一个通用的方法symplify,它试图用各种方法将表达式简化:

sympy.simplify(sympy.sin(x)**2+sympy.cos(x)**2)

[1 ]

sympy.simplify((x**3+x**2-x-1)/(x**2+2*x+1))

[x - 1 ]

sympy.simplify(sympy.gamma(x)/sympy.gamma(x-2))

[left(x - 2 ight) left(x - 1 ight) ]

simplify也有一个陷阱,因为它尝试用各种方法来找到最简化的结果,而“最简化的结果”本身就是一个不好定义的术语,比如,我们认为x** 2+2* xy+y** 2的最简化结果是(x+y)** 2

sympy.simplify(x**2+2*x*y+y**2)

[x^{2} + 2 x y + y^{2} ]

我们并没有得到我们想得到的(x+y)** 2的结果,针对这种情况,有一个方法factor,可以进行这样的简化。

sympy.factor(x**2+2*x*y+y**2)

[left(x + y ight)^{2} ]

simplify的第二个陷阱是它太慢了,因为它尝试各种方法来进行简化,如果我们事先知道最终简化的方式,最好选择特定的简化方法。使用simplify的最佳时机是交互式使用,来为进一步简化提供思路。

多项式/有理式函数简化

expand

expand,顾名思义,是展开:

sympy.expand((x+1)**3)

[x^{3} + 3 x^{2} + 3 x + 1 ]

sympy.expand((x+1)*(x+2)-(x-1)*x)

[4 x + 2 ]

factor

factor,是因式分解:

sympy.factor(x**2*z+4*x*y*z+4*y**2*z)

[z left(x + 2 y ight)^{2} ]

对于多项式,factorexpand的对立,factor可以保证每个因子都是不可约的,facor_list返回更有结构层次的输出:

sympy.factor_list(x**2+2*x*y+y**2)

[left ( 1, quad left [ left ( x + y, quad 2 ight ) ight ] ight ) ]

help(sympy.factor_list)
Help on function factor_list in module sympy.polys.polytools:

factor_list(f, *gens, **args)
    Compute a list of irreducible factors of ``f``.
    
    Examples
    ========
    
    >>> from sympy import factor_list
    >>> from sympy.abc import x, y
    
    >>> factor_list(2*x**5 + 2*x**4*y + 4*x**3 + 4*x**2*y + 2*x + 2*y)
    (2, [(x + y, 1), (x**2 + 1, 2)])
sympy.factor(x**2+2*x*y+y**2+y)

[x^{2} + 2 x y + y^{2} + y ]

sympy.factor_list(x**2+2*x*y*y**2+y)

[left ( 1, quad left [ left ( x^{2} + 2 x y^{3} + y, quad 1 ight ) ight ] ight ) ]

(2*x**5 + 2*x**4*y + 4*x**3 + 4*x**2*y + 2*x + 2*y)

[2 x^{5} + 2 x^{4} y + 4 x^{3} + 4 x^{2} y + 2 x + 2 y ]

sympy.factor(2*x**5 + 2*x**4*y + 4*x**3 + 4*x**2*y + 2*x + 2*y)

[2 left(x + y ight) left(x^{2} + 1 ight)^{2} ]

sympy.factor_list(2*x**5 + 2*x**4*y + 4*x**3 + 4*x**2*y + 2*x + 2*y)

[left ( 2, quad left [ left ( x + y, quad 1 ight ), quad left ( x^{2} + 1, quad 2 ight ) ight ] ight ) ]

collect

help(sympy.collect)
Help on function collect in module sympy.simplify.radsimp:

collect(expr, syms, func=None, evaluate=None, exact=False, distribute_order_term=True)
    Collect additive terms of an expression.
    
    This function collects additive terms of an expression with respect
    to a list of expression up to powers with rational exponents. By the
    term symbol here are meant arbitrary expressions, which can contain
    powers, products, sums etc. In other words symbol is a pattern which
    will be searched for in the expression's terms.
    
    The input expression is not expanded by :func:`collect`, so user is
    expected to provide an expression is an appropriate form. This makes
    :func:`collect` more predictable as there is no magic happening behind the
    scenes. However, it is important to note, that powers of products are
    converted to products of powers using the :func:`expand_power_base`
    function.
    
    There are two possible types of output. First, if ``evaluate`` flag is
    set, this function will return an expression with collected terms or
    else it will return a dictionary with expressions up to rational powers
    as keys and collected coefficients as values.
    
    Examples
    ========
    
    >>> from sympy import S, collect, expand, factor, Wild
    >>> from sympy.abc import a, b, c, x, y, z
    
    This function can collect symbolic coefficients in polynomials or
    rational expressions. It will manage to find all integer or rational
    powers of collection variable::
    
        >>> collect(a*x**2 + b*x**2 + a*x - b*x + c, x)
        c + x**2*(a + b) + x*(a - b)
    
    The same result can be achieved in dictionary form::
    
        >>> d = collect(a*x**2 + b*x**2 + a*x - b*x + c, x, evaluate=False)
        >>> d[x**2]
        a + b
        >>> d[x]
        a - b
        >>> d[S.One]
        c
    
    You can also work with multivariate polynomials. However, remember that
    this function is greedy so it will care only about a single symbol at time,
    in specification order::
    
        >>> collect(x**2 + y*x**2 + x*y + y + a*y, [x, y])
        x**2*(y + 1) + x*y + y*(a + 1)
    
    Also more complicated expressions can be used as patterns::
    
        >>> from sympy import sin, log
        >>> collect(a*sin(2*x) + b*sin(2*x), sin(2*x))
        (a + b)*sin(2*x)
    
        >>> collect(a*x*log(x) + b*(x*log(x)), x*log(x))
        x*(a + b)*log(x)
    
    You can use wildcards in the pattern::
    
        >>> w = Wild('w1')
        >>> collect(a*x**y - b*x**y, w**y)
        x**y*(a - b)
    
    It is also possible to work with symbolic powers, although it has more
    complicated behavior, because in this case power's base and symbolic part
    of the exponent are treated as a single symbol::
    
        >>> collect(a*x**c + b*x**c, x)
        a*x**c + b*x**c
        >>> collect(a*x**c + b*x**c, x**c)
        x**c*(a + b)
    
    However if you incorporate rationals to the exponents, then you will get
    well known behavior::
    
        >>> collect(a*x**(2*c) + b*x**(2*c), x**c)
        x**(2*c)*(a + b)
    
    Note also that all previously stated facts about :func:`collect` function
    apply to the exponential function, so you can get::
    
        >>> from sympy import exp
        >>> collect(a*exp(2*x) + b*exp(2*x), exp(x))
        (a + b)*exp(2*x)
    
    If you are interested only in collecting specific powers of some symbols
    then set ``exact`` flag in arguments::
    
        >>> collect(a*x**7 + b*x**7, x, exact=True)
        a*x**7 + b*x**7
        >>> collect(a*x**7 + b*x**7, x**7, exact=True)
        x**7*(a + b)
    
    You can also apply this function to differential equations, where
    derivatives of arbitrary order can be collected. Note that if you
    collect with respect to a function or a derivative of a function, all
    derivatives of that function will also be collected. Use
    ``exact=True`` to prevent this from happening::
    
        >>> from sympy import Derivative as D, collect, Function
        >>> f = Function('f') (x)
    
        >>> collect(a*D(f,x) + b*D(f,x), D(f,x))
        (a + b)*Derivative(f(x), x)
    
        >>> collect(a*D(D(f,x),x) + b*D(D(f,x),x), f)
        (a + b)*Derivative(f(x), x, x)
    
        >>> collect(a*D(D(f,x),x) + b*D(D(f,x),x), D(f,x), exact=True)
        a*Derivative(f(x), x, x) + b*Derivative(f(x), x, x)
    
        >>> collect(a*D(f,x) + b*D(f,x) + a*f + b*f, f)
        (a + b)*f(x) + (a + b)*Derivative(f(x), x)
    
    Or you can even match both derivative order and exponent at the same time::
    
        >>> collect(a*D(D(f,x),x)**2 + b*D(D(f,x),x)**2, D(f,x))
        (a + b)*Derivative(f(x), x, x)**2
    
    Finally, you can apply a function to each of the collected coefficients.
    For example you can factorize symbolic coefficients of polynomial::
    
        >>> f = expand((x + a + 1)**3)
    
        >>> collect(f, x, factor)
        x**3 + 3*x**2*(a + 1) + 3*x*(a + 1)**2 + (a + 1)**3
    
    .. note:: Arguments are expected to be in expanded form, so you might have
              to call :func:`expand` prior to calling this function.
    
    See Also
    ========
    collect_const, collect_sqrt, rcollect
expr=x*y+x-3+2*x**2-z*x**2+x**3;expr

[x^{3} - x^{2} z + 2 x^{2} + x y + x - 3 ]

sympy.collect(expr,x)

[x^{3} + x^{2} left(- z + 2 ight) + x left(y + 1 ight) - 3 ]

collect经常与coeff一起使用,来给出第n次幂的系数:

sympy.collect(expr,x).coeff(x,2)

[- z + 2 ]

cancle

通分及约分:

sympy.cancel((x**2+2*x+1)/(x**2+x))

[frac{1}{x} left(x + 1 ight) ]

expr=1/x+(3*x/2-2)/(x-4);expr

[frac{frac{3 x}{2} - 2}{x - 4} + frac{1}{x} ]

sympy.cancel(expr)

[frac{3 x^{2} - 2 x - 8}{2 x^{2} - 8 x} ]

apart

expr = (4*x**3 + 21*x**2 + 10*x + 12)/(x**4 + 5*x**3 + 5*x**2 + 4*x);expr

[frac{4 x^{3} + 21 x^{2} + 10 x + 12}{x^{4} + 5 x^{3} + 5 x^{2} + 4 x} ]

sympy.apart(expr)

[frac{2 x - 1}{x^{2} + x + 1} - frac{1}{x + 4} + frac{3}{x} ]

help(sympy.apart)
Help on function apart in module sympy.polys.partfrac:

apart(f, x=None, full=False, **options)
    Compute partial fraction decomposition of a rational function.
    
    Given a rational function ``f``, computes the partial fraction
    decomposition of ``f``. Two algorithms are available: One is based on the
    undertermined coefficients method, the other is Bronstein's full partial
    fraction decomposition algorithm.
    
    The undetermined coefficients method (selected by ``full=False``) uses
    polynomial factorization (and therefore accepts the same options as
    factor) for the denominator. Per default it works over the rational
    numbers, therefore decomposition of denominators with non-rational roots
    (e.g. irrational, complex roots) is not supported by default (see options
    of factor).
    
    Bronstein's algorithm can be selected by using ``full=True`` and allows a
    decomposition of denominators with non-rational roots. A human-readable
    result can be obtained via ``doit()`` (see examples below).
    
    Examples
    ========
    
    >>> from sympy.polys.partfrac import apart
    >>> from sympy.abc import x, y
    
    By default, using the undetermined coefficients method:
    
    >>> apart(y/(x + 2)/(x + 1), x)
    -y/(x + 2) + y/(x + 1)
    
    The undetermined coefficients method does not provide a result when the
    denominators roots are not rational:
    
    >>> apart(y/(x**2 + x + 1), x)
    y/(x**2 + x + 1)
    
    You can choose Bronstein's algorithm by setting ``full=True``:
    
    >>> apart(y/(x**2 + x + 1), x, full=True)
    RootSum(_w**2 + _w + 1, Lambda(_a, (-2*_a*y/3 - y/3)/(-_a + x)))
    
    Calling ``doit()`` yields a human-readable result:
    
    >>> apart(y/(x**2 + x + 1), x, full=True).doit()
    (-y/3 - 2*y*(-1/2 - sqrt(3)*I/2)/3)/(x + 1/2 + sqrt(3)*I/2) + (-y/3 -
        2*y*(-1/2 + sqrt(3)*I/2)/3)/(x + 1/2 - sqrt(3)*I/2)
    
    
    See Also
    ========
    
    apart_list, assemble_partfrac_list

三角几何简化

sympy与Python保持一致,在三角函数名前加a,就是相应的反三角函数

sympy.acos(x)

[operatorname{acos}{left (x ight )} ]

sympy.cos(sympy.acos(x))

[x ]

sympy.asin(1)

[frac{pi}{2} ]

trigsimp

sympy.trigsimp(sympy.sin(x)**2+sympy.cos(x)**2)

[1 ]

sympy.trigsimp(sympy.sin(x)**4-2*sympy.cos(x)**2*sympy.sin(x)**2+sympy.cos(x)**4)

[frac{1}{2} cos{left (4 x ight )} + frac{1}{2} ]

sympy.trigsimp(sympy.cosh(x)**2+sympy.sinh(x)**2)

[cosh{left (2 x ight )} ]

sympy.trigsimp(sympy.sinh(x)/sympy.cosh(x))

[ anh{left (x ight )} ]

expand_trig

sympy.expand_trig(sympy.sin(2*x))

[2 sin{left (x ight )} cos{left (x ight )} ]

sympy.expand_trig(sympy.tan(2*x))

[frac{2 an{left (x ight )}}{- an^{2}{left (x ight )} + 1} ]

expand_trig是试图把三角函数表达式变大,trigsimp是试图将三角函数表达式变小。

Powers

先看3个等式:

    1. xa*xb=x^(a+b)
  • 2.xa*ya=(xy)^a
  • 3.(xa)b=x^(ab)

等式1总是成立的;
等式2并不是总是成立的,比如,x=y=-1,a=1/2,左边为-1,右边为1
等式3也并不是总成立的,比如,x=-1,a=2,b=1/2,左边为-1,右边为1

如果幂表达式并不总是正确的,Sympy就不会进行简化,我们需要对符号进行做一些假定。

  • Sympy符号默认是复数,这意味着不会对一些给定的表达式进行简化,除非它适合所有的复数。
  • 符号可以通过传入symbols的参数来完成一些假定,对于本部分剩余的地方,我们假定x,y是正数,a,b是实数,z,t,c是任意的复数,来看看能得到什么结果。
x,y=sympy.symbols('x y',positive=True)
a,b=sympy.symbols('a b',real=True)
z,t,c=sympy.symbols('z t c')

需要注意的是,sqrt(x)x**Rational(1,2)的快捷方式:

sympy.sqrt(x)==x**sympy.Rational(1,2)
True

powsimp

该方法施加等式1,等式2

sympy.powsimp(x**a*x**b)

[x^{a + b} ]

sympy.powsimp(x**a*y**a)

[left(x y ight)^{a} ]

注意,powsimp不会对不合法的结果做简化,如

sympy.powsimp(t**c*z**c)

[t^{c} z^{c} ]

如果要强行进行简化,设置force=True:

sympy.powsimp(t**c*z**c,force=True)

[left(t z ight)^{c} ]

sympy.powsimp(z**2*t**2)

[t^{2} z^{2} ]

expand_power_exp/expand_power_base

sympy.expand_power_exp(x**(a+b)) # 简化指数

[x^{a} x^{b} ]

help(sympy.expand_power_exp)
Help on function expand_power_exp in module sympy.core.function:

expand_power_exp(expr, deep=True)
    Wrapper around expand that only uses the power_exp hint.
    
    See the expand docstring for more information.
    
    Examples
    ========
    
    >>> from sympy import expand_power_exp
    >>> from sympy.abc import x, y
    >>> expand_power_exp(x**(y + 2))
    x**2*x**y
sympy.expand_power_exp((x*y)**a) # 简化指数

[left(x y ight)^{a} ]

sympy.expand_power_base((x*y)**a) # 简化底

[x^{a} y^{a} ]

就像powsimp,一样,如果非法,expand_power_base不会进行简化:

sympy.expand_power_base((z*t)**c)

[left(t z ight)^{c} ]

但它也有force=True:

sympy.expand_power_base((z*t)**c,force=True)

[t^{c} z^{c} ]

powdenest

help(sympy.powdenest)
Help on function powdenest in module sympy.simplify.powsimp:

powdenest(eq, force=False, polar=False)
    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
sympy.powdenest((x**a)**b)

[x^{a b} ]

sympy.powdenest((z*a)**b)

[left(a z ight)^{b} ]

sympy
<module 'sympy' from 'E:\software\Anaconda\lib\site-packages\sympy\__init__.py'>

指数和对数

先看2个等式:

  • log(xy)=lob(x)+log(y)
  • log(x^n)=nlog(x)
    两个等式都不是总成立的,然而,如果x,y是正数,且n为实数,则等式成立。

在sympy中,log是自然对数,也就是ln,为了防止忘了,sympy自动提供log=ln.

sympy.log(x)

[log{left (x ight )} ]

sympy.ln(x)

[log{left (x ight )} ]

x,y=sympy.symbols('x y',positive=True)
n=sympy.symbols('n',real=True)

像之前的t,z,没有任何假定。

expand_log

expand_log来对两个等式,从左到右的计算:

sympy.expand_log(sympy.log(x*y))

[log{left (x ight )} + log{left (y ight )} ]

sympy.expand_log(sympy.log(x/y))

[log{left (x ight )} - log{left (y ight )} ]

sympy.log(x/y)

[log{left (frac{x}{y} ight )} ]

sympy.expand_log(sympy.log(x**2))

[2 log{left (x ight )} ]

sympy.expand_log(sympy.log(x**n))

[n log{left (x ight )} ]

sympy.expand_log(sympy.log(t*z))

[log{left (t z ight )} ]

当然,expand_log也有个force=True设置,可以强制的进行log 展开:

sympy.expand_log(sympy.log(t*z),force=True)

[log{left (t ight )} + log{left (z ight )} ]

logcombine

对这两个等式,从右往左计算:

sympy.logcombine(sympy.log(x)+sympy.log(y))

[log{left (x y ight )} ]

sympy.logcombine(sympy.log(x)-sympy.log(y))

[log{left (frac{x}{y} ight )} ]

同样,logcombine也有force=True的设置

特殊函数

sympy有很多个特殊函数,从组合数学到数学物理方程,详尽的描述见https://docs.sympy.org/latest/modules/functions/index.html#functions-contents

首先定义x,y,z,k,m,l,对它们不做任何假定。

x,y,z,k,m,l=sympy.symbols('x y z k m l')

n的阶乘用fatorial(n):

sympy.factorial(n)

[n! ]

二项式系数方程,binominal(n,k)

sympy.binomial(n,k)

[{inom{n}{k}} ]

gamma函数,gamma(z)

sympy.gamma(z)

[Gamma{left(z ight)} ]

重写rewrite

对于特殊函数的常规操作是重写,这实际上对于sympy中所有函数都成立,并不仅仅局限于特殊函数,用expr.rewrite(function)来完成重写。

sympy.tan(x).rewrite(sympy.sin)

[frac{2 sin^{2}{left (x ight )}}{sin{left (2 x ight )}} ]

sympy.factorial(x).rewrite(sympy.gamma)

[Gamma{left(x + 1 ight)} ]

expand_func

sympy.expand_func(sympy.gamma(3+x))

[x left(x + 1 ight) left(x + 2 ight) Gamma{left(x ight)} ]

combsimp

n,k=sympy.symbols('n k',integer=True)
sympy.combsimp(sympy.factorial(n)/sympy.factorial(n-3))

[n left(n - 2 ight) left(n - 1 ight) ]

sympy.combsimp(sympy.binomial(n+1,k+1)/sympy.binomial(n,k))

[frac{n + 1}{k + 1} ]

例子(连续分式)

def list_to_fac(l):
    expr=sympy.Integer(0)
    for i in reversed(l[1:]):
        expr+=i
        expr=1/expr
    return l[0]+expr
        
list_to_fac([x,y,z])

[x + frac{1}{y + frac{1}{z}} ]

我们在list_to_fac中用了sympy.Integer(0),所以就算传入list_to_fac的列表中都是Python int,最后返回结果还是sympy 对象。

list_to_fac([1,2,3,4])

[frac{43}{30} ]

sysm=sympy.symbols('a0:5')
a0, a1, a2, a3, a4 = sysm
fac=list_to_fac(sysm);fac

[a_{0} + frac{1}{a_{1} + frac{1}{a_{2} + frac{1}{a_{3} + frac{1}{a_{4}}}}} ]

sympy.cancel(fac)

[frac{a_{0} a_{1} a_{2} a_{3} a_{4} + a_{0} a_{1} a_{2} + a_{0} a_{1} a_{4} + a_{0} a_{3} a_{4} + a_{0} + a_{2} a_{3} a_{4} + a_{2} + a_{4}}{a_{1} a_{2} a_{3} a_{4} + a_{1} a_{2} + a_{1} a_{4} + a_{3} a_{4} + 1} ]

fac=sympy.cancel(fac)
fac

[frac{a_{0} a_{1} a_{2} a_{3} a_{4} + a_{0} a_{1} a_{2} + a_{0} a_{1} a_{4} + a_{0} a_{3} a_{4} + a_{0} + a_{2} a_{3} a_{4} + a_{2} + a_{4}}{a_{1} a_{2} a_{3} a_{4} + a_{1} a_{2} + a_{1} a_{4} + a_{3} a_{4} + 1} ]

l=[]
fac=sympy.apart(fac,a1);fac

[a_{0} + frac{a_{2} a_{3} a_{4} + a_{2} + a_{4}}{a_{1} a_{2} a_{3} a_{4} + a_{1} a_{2} + a_{1} a_{4} + a_{3} a_{4} + 1} ]

l.append(a0)
fac = 1/(fac - a0);fac

[frac{1}{a_{2} a_{3} a_{4} + a_{2} + a_{4}} left(a_{1} a_{2} a_{3} a_{4} + a_{1} a_{2} + a_{1} a_{4} + a_{3} a_{4} + 1 ight) ]

fac=sympy.apart(fac,a1);fac

[a_{1} + frac{a_{3} a_{4} + 1}{a_{2} a_{3} a_{4} + a_{2} + a_{4}} ]

l.append(a1)

微积分

导数

sympy.diff(x**2,x)

[2 x ]

diff可以一次性求解多阶导数,比如求解n阶导数,可以将n个x传入diff中或者将x,n传入diff中:

sympy.diff(sympy.exp(x**2),x,2)

[2 left(2 x^{2} + 1 ight) e^{x^{2}} ]

sympy.diff(sympy.exp(x**2),x,x,)

[2 left(2 x^{2} + 1 ight) e^{x^{2}} ]

diff也可以同时对多个变量求导,将需要求导的变量及次数依次输入:

expr=sympy.exp(x*y*z)
expr

[e^{x y z} ]

expr.diff(x,y,y,z,z,z,z)

[x^{3} y^{2} left(x^{3} y^{3} z^{3} + 14 x^{2} y^{2} z^{2} + 52 x y z + 48 ight) e^{x y z} ]

用类Derivative来表示没有计算的导数,偏导数,与diff有相同的语法:

sympy.Derivative(expr,x,y,2,z,4)

[frac{partial^{7}}{partial xpartial y^{2}partial z^{4}} e^{x y z} ]

计算用类Derivative表示的导数,偏导数,用方法doit():

sympy.Derivative(expr,x,y,2,z,4).doit()

[x^{3} y^{2} left(x^{3} y^{3} z^{3} + 14 x^{2} y^{2} z^{2} + 52 x y z + 48 ight) e^{x y z} ]

不确定是几阶导数的可以用(x,n)元组传入diff:

m,n,a,b=sympy.symbols('m n a b')
expr=(a*x+b)**m;expr

[left(a x + b ight)^{m} ]

# sympy.diff(expr,(x,n))  # go wrong?!

积分

integrate来进行积分计算,包括定积分,不定积分,计算不定积分,只要把自变量传入进去就可以了。

sympy.integrate(sympy.cos(x),x)

[sin{left (x ight )} ]

如果要求定积分,传入(integration_variable, lower_limit, upper_limit):

from sympy import *
sympy.integrate(sympy.exp(-x),(x,0,oo))

[1 ]

sympy.integrate(sympy.exp(-x**2-y**2),(x,-oo,+oo),(y,-oo,oo))

[pi ]

如果integrate不能计算出一个积分值,会返回一个没有被计算的Integrate对象。

expr=sympy.integrate(x**x,x)
expr

[int x^{x}\, dx ]

Derivative一样,Integraral返回一个没被计算的积分对象,然后可以用doit来计算。

expr=sympy.Integral(sympy.log(x)**2,x)
expr

[int log^{2}{left (x ight )}\, dx ]

expr.doit()

[x log^{2}{left (x ight )} - 2 x log{left (x ight )} + 2 x ]

integ = Integral((x**4 + x**2*exp(x) - x**2 - 2*x*exp(x) - 2*x -
    exp(x))*exp(x)/((x - 1)**2*(x + 1)**2*(exp(x) + 1)), x)
integ

[int frac{left(x^{4} + x^{2} e^{x} - x^{2} - 2 x e^{x} - 2 x - e^{x} ight) e^{x}}{left(x - 1 ight)^{2} left(x + 1 ight)^{2} left(e^{x} + 1 ight)}\, dx ]

integ.doit()

[log{left (e^{x} + 1 ight )} + frac{e^{x}}{x^{2} - 1} ]

integ = Integral(sin(x**2), x);integ

[int sin{left (x^{2} ight )}\, dx ]

integ.doit()

[frac{3 sqrt{2} sqrt{pi} Sleft(frac{sqrt{2} x}{sqrt{pi}} ight)}{8 Gamma{left(frac{7}{4} ight)}} Gamma{left(frac{3}{4} ight)} ]

integ = Integral(x**y*exp(-x), (x, 0, oo));integ

[int_{0}^{infty} x^{y} e^{- x}\, dx ]

integ.doit()

[egin{cases} Gamma{left(y + 1 ight)} & ext{for}: - Re{y} < 1 \int_{0}^{infty} x^{y} e^{- x}\, dx & ext{otherwise} end{cases} ]

最后一个返回了一个Piecewise(分段函数)对象

求极限

求极限的基本语法:
limit(f(x), x, x0)

limit(sin(x)/x, x, 0)

[1 ]

当表达式在某点出奇异点时候,应该用limit,而不是subs

expr=x**2/exp(x);expr

[x^{2} e^{- x} ]

expr.subs(x,oo)

[mathrm{NaN} ]

limit(expr,x,oo)

[0 ]

limit(expr,x,-oo)

[infty ]

极限也有自己的对象,Limit,可以用doit来对该对象进行计算:

expr=Limit((cos(x)-1)/x,x,0);expr

[lim_{x o 0^+}left(frac{1}{x} left(cos{left (x ight )} - 1 ight) ight) ]

expr.doit()

[0 ]

如果仅仅要计算某一边的极限,最后再传入'+',或者'-'即可:

limit(1/x,x,0,'+')

[infty ]

limit(1/x,x,0,'-')

[-infty ]

级数展开

sympy可以在某点对一个函数进行非对称级数展开。如,为了计算f(x)在x=x0处,x^n阶的展开,用f(x).series(x,x0,n),x0和n可以省略,这种情况,默认x0=0,n=6.

expr=exp(sin(x));expr

[e^{sin{left (x ight )}} ]

expr.series(x,0,4)

[1 + x + frac{x^{2}}{2} + mathcal{O}left(x^{4} ight) ]

如果不保留高阶项,用removeO(Big O):

expr.series(x,0,4).removeO()

[frac{x^{2}}{2} + x + 1 ]

有限差分

differentiate_finite来求有限差分:

f,g=sympy.symbols('f g',cls=Function)
differentiate_finite(f(x)*g(x))
---------------------------------------------------------------------------

NameError                                 Traceback (most recent call last)

<ipython-input-210-faa797f23aa1> in <module>()
----> 1 differentiate_finite(f(x)*g(x))


NameError: name 'differentiate_finite' is not defined

设置参数evaluate=True来展开中间导数:

differentiate_finite(f(x)*g(x),evaluate=True)

如果已经有了Derivative对象,可以用as_finite_difference方法,来得到任意阶的近似的导数。

f=Function('f')
dfdx=f(x).diff(x)
dfdx
dfdx.as_finite_difference()

这里,用最少的点(对于1阶导数是2个点),在x处近似的得到f(x)的导数;

f(x)
dfdx=f(x).diff(x,2);dfdx
h=symbols('h')
pp=dfdx.as_finite_difference([-3*h,h,2*h]);pp
help(dfdx.as_finite_difference)
finite_diff_weights(2, [-3, -1, 2], 0)[-1][-1]

解方程

在sympy中,等号用Eq表示:

Eq(x,y)

然而,还有更简单的方法表示等式,在sympy中,凡是不是用Eq的表达式,自动的都认为是等于0:

solveset(Eq(x**2,1),x)
solveset(Eq(x**2-1),x)
solveset(Eq(x**2-1,0),x)
solveset(x**2-1,x)

解代数方程

解代数方程的主要方法是solveset,语法为:solveset(equation, variable=None, domain=S.Complexes),其中equation可以是Eq对象,或者默认为0的表达式。
注意还有个solve方程,语法为:solve(equations, variables) ,还是建议先用solveset

solveset返回值是 FiniteSet , Interval 或者 ImageSet 对象:

solveset(x**2 - x, x)
solveset(x - x, x, domain=S.Reals)
solveset(sin(x) - 1, x, domain=S.Reals)
solveset(exp(x), x) 
solveset(cos(x) - x, x) 

在代数方程中,如果是线性方程,现在用linsolve,以后的版本可能会直接用solveset:

help(linsolve)
  • 方程形式的列表:
linsolve([x + y + z - 1, x + y + 2*z - 3 ], (x, y, z))
  • 增广矩阵形式:
linsolve(Matrix(([1, 1, 1, 1], [1, 1, 2, 3])), (x, y, z))
  • A * x = b形式
M = Matrix(((1, 1, 1, 1), (1, 1, 2, 3)))
M
system = A, b = M[:, :-1], M[:, -1]
system
linsolve(system, x, y, z)

在模块solveset中,用nonlinsolve来求解非线性方程组:

  • 仅有实数解
a, b, c, d = symbols('a, b, c, d', real=True)
nonlinsolve([a**2 + a, a - b], [a, b]) # that means,a**2+a=0;a-b=0
  • 仅有复数解
nonlinsolve([x**2 + 1, y**2 + 1], [x, y])
help(nonlinsolve)
  • 既有实根也有复数根
system = [x**2 - 2*y**2 -2, x*y - 2]
vars = [x, y]
nonlinsolve(system, vars)
  • 无限个解
nonlinsolve([x*y, x*y - x], [x, y])

注意
(1):解的顺序对应了symbols的顺序
(2)当前,nonlinsolve不返回以LambertW形式的解

solve([x**2 - y**2/exp(x)], [x, y], dict=True)
solve([sin(x + y), cos(x - y)], [x, y]) #可用于这种情况,但不会给出全部解

solveset只返回零点一次,用roots可以得到零点及相应出现的次数。

solveset(x**2,x)
roots(x**2,x)

解微分方程

dsolve解微分方程。

f, g = symbols('f g', cls=Function)
f(x).diff(x,2)
diffeq=Eq(f(x).diff(x,2)-2*f(x).diff(x)+f(x),sin(x));diffeq
dsolve(diffeq,f(x))

小总结
关于各种solve的传入参数,第一个是方程(组),第二个是相应的未知量(有必要写出来,它的顺序决定了求解值的顺序)

矩阵

Matrix对象来表示矩阵。

M=Matrix([[1,2,3],[2,3,4]]);M

为了简便,只传入一个列表,M即代表是列向量。

Matrix([1,2,3])
N = Matrix([0, 1, 1])
M*N

关于sympy中的矩阵一个比较重要的点是:不像sympy中其他对象,Matrix是可变的。这意味着,它们可以原地改变。这样的一个缺点是,它不能用在需要是不可变对象的地方,比如在sympy表达式中,如果需要一个Matrix的不变的版本,用ImmutableMatrix.

基本操作

  • shape
M
M.shape
  • 查询行列
M.row(0)
M.col(1)
  • 增删行列

删除一行或者一列,用row_del,column_del,这些操作都发生在原地。

M

[left[egin{matrix}3 & -2 & 4 & -2\5 & 3 & -3 & -2\5 & -2 & 2 & -2\5 & -2 & -3 & 3end{matrix} ight] ]

M.col_del(0)
M

[left[egin{matrix}-2 & 4 & -2\3 & -3 & -2\-2 & 2 & -2\-2 & -3 & 3end{matrix} ight] ]

M.row_del(1)
M

[left[egin{matrix}-2 & 4 & -2\-2 & 2 & -2\-2 & -3 & 3end{matrix} ight] ]

增加行列,用row_insert,col_insert,这些操作不会在原地进行。

M.row_insert(1,Matrix([[0,1]]))
---------------------------------------------------------------------------

ShapeError                                Traceback (most recent call last)

<ipython-input-270-8ea1d073fb8d> in <module>()
----> 1 M.row_insert(1,Matrix([[0,1]]))


E:softwareAnacondalibsite-packagessympymatricesmatrices.py in row_insert(self, pos, mti)
   4039         if self.cols != mti.cols:
   4040             raise ShapeError(
-> 4041                 "`self` and `mti` must have the same number of columns.")
   4042 
   4043         newmat = self.zeros(self.rows + mti.rows, self.cols)


ShapeError: `self` and `mti` must have the same number of columns.
M

[left[egin{matrix}1 & 3\0 & 1\-2 & 3end{matrix} ight] ]

M=M.row_insert(1,Matrix([[0,1]]))
M = M.col_insert(0, Matrix([1, -2])) # 这里,列向量,只需要一个list,不用嵌套。
---------------------------------------------------------------------------

ShapeError                                Traceback (most recent call last)

<ipython-input-280-3d900b7dc745> in <module>()
----> 1 M = M.col_insert(0, Matrix([1, -2])) # 这里,列向量,只需要一个list,不用嵌套。


E:softwareAnacondalibsite-packagessympymatricesmatrices.py in col_insert(self, pos, mti)
   4070         """
   4071         if pos == 0:
-> 4072             return mti.row_join(self)
   4073         elif pos < 0:
   4074             pos = self.cols + pos


E:softwareAnacondalibsite-packagessympymatricesmatrices.py in row_join(self, rhs)
   3966         if self.rows != rhs.rows:
   3967             raise ShapeError(
-> 3968                 "`self` and `rhs` must have the same number of rows.")
   3969 
   3970         from sympy.matrices import MutableMatrix


ShapeError: `self` and `rhs` must have the same number of rows.
M

[left[egin{matrix}-2 & 4 & -2\-2 & 2 & -2\-2 & -3 & 3end{matrix} ight] ]

除非有特殊声明,否则下面的操作都不会在原地进行,通常,不在原地进行的操作会返回一个新的对象,在原地操作会返回None.

基本操作

M = Matrix([[1, 3], [-2, 3]])
N = Matrix([[0, 3], [0, 7]])
M

[left[egin{matrix}1 & 3\-2 & 3end{matrix} ight] ]

N

[left[egin{matrix}0 & 3\0 & 7end{matrix} ight] ]

M+N

[left[egin{matrix}1 & 6\-2 & 10end{matrix} ight] ]

M*2

[left[egin{matrix}2 & 6\-4 & 6end{matrix} ight] ]

M**2

[left[egin{matrix}-5 & 12\-8 & 3end{matrix} ight] ]

M*M  # 与上式一样

[left[egin{matrix}-5 & 12\-8 & 3end{matrix} ight] ]

M**-1

[left[egin{matrix}frac{1}{3} & - frac{1}{3}\frac{2}{9} & frac{1}{9}end{matrix} ight] ]

M**-1*M

[left[egin{matrix}1 & 0\0 & 1end{matrix} ight] ]

M.T

[left[egin{matrix}1 & -2\3 & 3end{matrix} ight] ]

Matrix constructor

eye(3)

[left[egin{matrix}1 & 0 & 0\0 & 1 & 0\0 & 0 & 1end{matrix} ight] ]

eye(4)

[left[egin{matrix}1 & 0 & 0 & 0\0 & 1 & 0 & 0\0 & 0 & 1 & 0\0 & 0 & 0 & 1end{matrix} ight] ]

zeros(2,3)

[left[egin{matrix}0 & 0 & 0\0 & 0 & 0end{matrix} ight] ]

ones(2,3)

[left[egin{matrix}1 & 1 & 1\1 & 1 & 1end{matrix} ight] ]

对角阵用diag,其参数既可以是数字也可以是矩阵,剩下的是0.

diag(1,2,3)

[left[egin{matrix}1 & 0 & 0\0 & 2 & 0\0 & 0 & 3end{matrix} ight] ]

diag(-1,ones(2,3),Matrix([2,3,4]))

[left[egin{matrix}-1 & 0 & 0 & 0 & 0\0 & 1 & 1 & 1 & 0\0 & 1 & 1 & 1 & 0\0 & 0 & 0 & 0 & 2\0 & 0 & 0 & 0 & 3\0 & 0 & 0 & 0 & 4end{matrix} ight] ]

高级方法

  • 行列式,用det
M = Matrix([[1, 0, 1], [2, -1, 3], [4, 3, 2]])
M

[left[egin{matrix}1 & 0 & 1\2 & -1 & 3\4 & 3 & 2end{matrix} ight] ]

M.det()

[-1 ]

  • RREF
help(M.rref)
Help on method rref in module sympy.matrices.matrices:

rref(iszerofunc=<function _iszero at 0x000001AAC40281E0>, simplify=False) method of sympy.matrices.dense.MutableDenseMatrix instance
    Return reduced row-echelon form of matrix and indices of pivot vars.
    
    To simplify elements before finding nonzero pivots set simplify=True
    (to use the default SymPy simplify function) or pass a custom
    simplify function.
    
    Examples
    ========
    
    >>> from sympy import Matrix
    >>> from sympy.abc import x
    >>> m = Matrix([[1, 2], [x, 1 - 1/x]])
    >>> m.rref()
    (Matrix([
    [1, 0],
    [0, 1]]), [0, 1])

rref 就是reduced row echelon form,即简化列梯形矩阵,返回值第一个是经过线性变换后的矩阵,第二个元组,包含了对角的列编号。

M = Matrix([[1, 0, 1, 3], [2, 3, 4, 7], [-1, -3, -3, -4]])
M

[left[egin{matrix}1 & 0 & 1 & 3\2 & 3 & 4 & 7\-1 & -3 & -3 & -4end{matrix} ight] ]

M.rref()

[left ( left[egin{matrix}1 & 0 & 1 & 3\0 & 1 & frac{2}{3} & frac{1}{3}\0 & 0 & 0 & 0end{matrix} ight], quad left [ 0, quad 1 ight ] ight ) ]

  • Nullspace

对于线性方程:Ax=b。如果有唯一解的化,那么没有任何向量乘以A后能够等于0向量:Ab=0。这个时候叫做full rank。也就是rank =n。
如果如果Ax=b没有唯一解,那么就会有很多向量使Ab=0。这些向量b本身组成了一个空间(可以哦嗯一组基表示所有满足要求的向量)。
如果rank=n-1,那么这个nullspace由一个基向量就可以表示,如果rank=n-2,那么这个nullspace由两个基向量才能表示。
所以求一个矩阵的nullspace也就是求这个space的基。
如果A是方阵,那么A的rank就等于不为零的特征值的个数。

M = Matrix([[1, 2, 3, 0, 0], [4, 10, 0, 0, 1]])
M

[left[egin{matrix}1 & 2 & 3 & 0 & 0\4 & 10 & 0 & 0 & 1end{matrix} ight] ]

M.rref()

[left ( left[egin{matrix}1 & 0 & 15 & 0 & -1\0 & 1 & -6 & 0 & frac{1}{2}end{matrix} ight], quad left [ 0, quad 1 ight ] ight ) ]

M.nullspace()

[left [ left[egin{matrix}-15\6\1\0\0end{matrix} ight], quad left[egin{matrix}0\0\0\1\0end{matrix} ight], quad left[egin{matrix}1\- frac{1}{2}\0\0\1end{matrix} ight] ight ] ]

Columnspace

M = Matrix([[1, 1, 2], [2 ,1 , 3], [3 , 1, 4]])

特征值,特征向量,对角化

M = Matrix([[3, -2,  4, -2], [5,  3, -3, -2], [5, -2,  2, -2], [5, -2, -3,  3]])
M

[left[egin{matrix}3 & -2 & 4 & -2\5 & 3 & -3 & -2\5 & -2 & 2 & -2\5 & -2 & -3 & 3end{matrix} ight] ]

M.eigenvals()

[left { -2 : 1, quad 3 : 1, quad 5 : 2 ight } ]

M.eigenvects()

[left [ left ( -2, quad 1, quad left [ left[egin{matrix}0\1\1\1end{matrix} ight] ight ] ight ), quad left ( 3, quad 1, quad left [ left[egin{matrix}1\1\1\1end{matrix} ight] ight ] ight ), quad left ( 5, quad 2, quad left [ left[egin{matrix}1\1\1\0end{matrix} ight], quad left[egin{matrix}0\-1\0\1end{matrix} ight] ight ] ight ) ight ] ]

P, D = M.diagonalize()
P

[left[egin{matrix}0 & 1 & 1 & 0\1 & 1 & 1 & -1\1 & 1 & 1 & 0\1 & 1 & 0 & 1end{matrix} ight] ]

D

[left[egin{matrix}-2 & 0 & 0 & 0\0 & 3 & 0 & 0\0 & 0 & 5 & 0\0 & 0 & 0 & 5end{matrix} ight] ]

P*D*P**-1

[left[egin{matrix}3 & -2 & 4 & -2\5 & 3 & -3 & -2\5 & -2 & 2 & -2\5 & -2 & -3 & 3end{matrix} ight] ]

如果想要的特征多项式,直接用charpoly,这比engenvals更高效,有时候符号表示的根求解非常麻烦。

lamda = symbols('lamda')  # lambda 是python中的关键字,所以sympy中用lamda
p = M.charpoly(lamda)
p

[operatorname{PurePoly}{left( lambda^{4} - 11 lambda^{3} + 29 lambda^{2} + 35 lambda - 150, lambda, domain=mathbb{Z} ight)} ]

factor(p)

[left(lambda - 5 ight)^{2} left(lambda - 3 ight) left(lambda + 2 ight) ]

高级表达式操作

理解表达式树

首先来看表达式在sympy中是怎么表达的,比如2^x+xy,可以用srepr来内在的看表达式是什么样子。

expr=2**x+x*y
expr

[2^{x} + x y ]

srepr(expr)
"Add(Pow(Integer(2), Symbol('x')), Mul(Symbol('x'), Symbol('y')))"

免责声明:文章转载自《sympy简明用法》仅用于学习参考。如对内容有疑问,请及时联系本站处理。

上篇Underlying cause: java.sql.SQLException : Access denied for user 'root'@'s150' (using password: YES)服务管理程序,检测服务是否存在,如果不存在,启动它下篇

宿迁高防,2C2G15M,22元/月;香港BGP,2C5G5M,25元/月 雨云优惠码:MjYwNzM=

相关文章

【Ray Tracing The Next Week 超详解】 光线追踪2-7 任意长方体 &amp;amp;&amp;amp; 场景案例

上一篇比较简单,很久才发是因为做了一些好玩的场景,后来发现这一章是专门写场景例子的,所以就安排到了这一篇 Preface 这一篇要介绍的内容有: 1. 自己做的光照例子 2. Cornell box画质问题及优化方案 3. 新的场景几何体——长方体 轴平行长方体 任意长方体 我们这一篇重实践轻理论阐述 ready 1. 需要上一章的知识 但是,上一章的Co...

全面总结:matlab怎么做漂亮的图

源地址:http://blog.csdn.net/ccxcau/article/details/7362764 MATLAB受到控制界广泛接受的一个重要原因是因为它提供了方便的绘图功能.本章主要介绍2维图形对象的生成函数及图形控制函数的使用方法,还将简单地介绍一 些图形的修饰与标注函数及操作和控制MATLAB各种图形对象的方法. 第一节 图形窗口与坐标系...

matlab绘图

MatLab绘图 作为一个功能强大的工具软件,Matlab具有很强的图形处理功能,提供了大量的二维、三维图形函数。由于系统采用面向对象的技术和丰富的矩阵运算,所以在图形处理方面方便又高效。 一般来说,一个命令行输入一条命令,命令行以回车结束。但一个命令行也可以输入若干条命令,各命令之间以逗号分隔,若前一命令后带有分号,则逗号可以省略。 如果一个命令行很长...

(笔记)Linux下C语言实现静态IP地址,掩码,网关的设置

#include <sys/ioctl.h> #include <sys/types.h> #include <sys/socket.h> #include <netinet/in.h> #include <arpa/inet.h> #include <net/if.h> #inclu...

matlab绘图方法汇总

Matlab画图 强大的画图功能是Matlab的特点之中的一个,Matlab提供了一系列的画图函数,用户不须要过多的考虑画图的细节,仅仅须要给出一些基本參数就能得到所需图形,这类函数称为高层画图函数。此外,Matlab还提供了直接对图形句柄进行操作的低层画图操作。这类操作将图形的每个图形元素(如坐标轴、曲线、文字等)看做一个独立的对象,系统给每个对象分配一...

深度学习入门|第五章 误差反向传播法

误差反向传播法 前言 此为本人学习《深度学习入门》的学习笔记,详情请阅读原书 数值微分虽然简单,也容易实现,但是缺点是计算上比较费时间,本章介绍一个高效计算权重参数的梯度的方法--误差反向传播法 一、计算图 计算图将计算过程用图形表示出来。这里说的图形是数据结构图,通过多个节点和边表示(连接节点的直线称为“边”)。 1、用计算图求解 实例:太郎在超市买了...