SymPy是一个用来处理数学符号的Python库,旨在成为一个多功能但代码尽可能简洁以便于理解和扩展的计算机代数系统(CAS)。同时,SymPy完全是用Python编写的,并且不依赖任何外部的库。
本教程主要讲述SymPy的基础知识,阅读它有助于你理解能用SymPy来做什么(怎样做),如果你想要更深入地学习SymPy,请阅读 SymPy User’s Guide, SymPy Modules Reference, 或者直接阅读 源代码 。
下载SymPy最简单的方法就是从 http://code.google.com/p/sympy/ 下载最新的压缩包:
解压:
$ tar xzf sympy-0.5.12.tar.gz
然后在Python解释器中尝试使用:
$ cd sympy-0.5.12 $ python Python 2.4.4 (#2, Jan 3 2008, 13:36:28) [GCC 4.2.3 20071123 (prerelease) (Debian 4.2.2-4)] on linux2 Type "help", "copyright", "credits" or "license" for more information. >>> from sympy import Symbol, cos >>> x = Symbol("x") >>> (1/cos(x)).series(x, 0, 10) 1 + x**2/2 + 5*x**4/24 + 61*x**6/720 + 277*x**8/8064 + O(x**10)
我们推荐你在你的程序中按上面的方法使用SymPy。你还可以像安装其它Python模块一样,通过 ./setup.py install 安装SymPy,或者就像安装一个软件包一样,在你喜欢的Linux发行版中安装它:
Installing SymPy in Debian
$ sudo apt-get install python-sympy Reading package lists... Done Building dependency tree Reading state information... Done The following NEW packages will be installed: python-sympy 0 upgraded, 1 newly installed, 0 to remove and 18 not upgraded. Need to get 991kB of archives. After this operation, 5976kB of additional disk space will be used. Get:1 http://ftp.cz.debian.org unstable/main python-sympy 0.5.12-1 [991kB] Fetched 991kB in 2s (361kB/s) Selecting previously deselected package python-sympy. (Reading database ... 232619 files and directories currently installed.) Unpacking python-sympy (from .../python-sympy_0.5.12-1_all.deb) ... Setting up python-sympy (0.5.12-1) ...
想了解其它安装SymPy的方法,参考wiki页面 Download and Installation 。
为了验证新特性,或者想弄清楚怎样使用这些特性,你可以使用一个对IPython封装好的包 isympy (如果从源代码目录中运行,则为 bin/isympy ),这只是一个已经导入了相关SymPy模块并且定义了符号x,y,z和一些其他东西的标准Python脚本。
$ cd sympy $ ./bin/isympy IPython console for SymPy 0.7.2-git (Python 2.7.1) (ground types: gmpy) These commands were executed: >>> from __future__ import division >>> from sympy import * >>> x, y, z, t = symbols('x y z t') >>> k, m, n = symbols('k m n', integer=True) >>> f, g, h = symbols('f g h', cls=Function) Documentation can be found at http://www.sympy.org In [1]: (1/cos(x)).series(x, 0, 10) Out[1]: 2 4 6 8 x 5*x 61*x 277*x / 10\ 1 + -- + ---- + ----- + ------ + O\x / 2 24 720 8064
注解
你所输入的命令是很粗糙的,所以我们在常规的Python解释器中要用三条命令才能完成的事,在isympy中用一条语句就可以完成。
SymPy有三个内置的数值类型:浮点型,分数型和整型
分数类通过一对整型数来表示一个分数:分子和分母,所以Rational(1, 2)表示1/2,Rational(5, 2)表示5/2,等等。
>>> from sympy import Rational
>>> a = Rational(1, 2)
>>> a
1/2
>>> a*2
1
>>> Rational(2)**50/Rational(10)**50
1/88817841970012523233890533447265625
因为你可能建立的是一个基于Python数值类型的数据,而不是基于SymPy数值类型的,所以请小心使用Python中的整型数和浮点数,特别是在除法运算中。两个Python中的整型数的比值可能会形成一个浮点数——这是Python 3所使用的”真正的除法”的标准,也是从__future__中导入除法的 isympy 的默认操作。
>>> 1/2
0.5
但是在一些早期的Python版本当中,除法还没有被导入,所以会对除法操作的结果取整,截断小数点后面的部分。
>>> 1/2
0
然而,不管在哪种情况下,你都不是在处理一个SymPy的数据,因为是Python创建它自己的数据。在大多数情况下,你可能需要处理分数,所以这时要确保SymPy的输出结果也是分数。当 R 等同于 Rational 时,更便于我们编写语句:
>>> R = Rational
>>> R(1, 2)
1/2
>>> R(1)/2 # R(1) is a SymPy Integer and Integer/int gives a Rational
1/2
SymPy中还有一些特殊的常量,例如e和pi,这些常量在SymPy中被看作是符号(1 + pi不会等于一个数值,它只是等于1 + pi),并且拥有任意精度:
>>> from sympy import pi, E
>>> pi**2
pi**2
>>> pi.evalf()
3.14159265358979
>>> (pi + E).evalf()
5.85987448204884
正如你所看到的,通过evalf函数对以上表达式求值,可以得到一个浮点数。
符号 oo 表示数学上所说的无穷:
>>> from sympy import oo
>>> oo > 99999
True
>>> oo + 1
oo
跟其它计算机代数系统不同,在SymPy中必须明确声明符号变量:
>>> from sympy import Symbol
>>> x = Symbol('x')
>>> y = Symbol('y')
左边的是常规的Python变量,已经将SymPy符号类分配给这些变量。预定义符号(包括含有希腊字母的符号)可以从abc中导入:
>>> from sympy.abc import x, theta
符号还可以通过 symbols 或者 var 函数来创建,其中,使用 var 函数会自动添加新建的符号到命名空间,这两个函数都允许使用范围表示作为参数:
>>> from sympy import symbols, var
>>> a, b, c = symbols('a,b,c')
>>> d, e, f = symbols('d:f')
>>> var('g:h')
(g, h)
>>> var('g:2')
(g0, g1)
符号类的实例能够”很好地结合”,由多个符号可以构建一个表达式:
>>> x + y + x - y
2*x
>>> (x + y)**2
(x + y)**2
>>> ((x + y)**2).expand()
x**2 + 2*x*y + y**2
使用函数 subs(old, new) 可以将给定的符号代换成其它的数字,符号或者表达式:
>>> ((x + y)**2).subs(x, 1)
(y + 1)**2
>>> ((x + y)**2).subs(x, y)
4*y**2
>>> ((x + y)**2).subs(x, 1 - y)
1
对于此教程的余下部分,我们假设已经运行以下语句:
>>> from sympy import init_printing
>>> init_printing(use_unicode=False, wrap_line=False, no_global=True)
运行这些语句会让输出的结果更雅观,查看以下章节 打印 ,如果你已经安装了unicode字体,可以设置 use_unicode=True,使得输出更加美观。
实现部分分数分解,使用函数 apart(expr, x) :
>>> from sympy import apart
>>> from sympy.abc import x, y, z
>>> 1/( (x + 2)*(x + 1) )
1
---------------
(x + 1)*(x + 2)
>>> apart(1/( (x + 2)*(x + 1) ), x)
1 1
- ----- + -----
x + 2 x + 1
>>> (x + 1)/(x - 1)
x + 1
-----
x - 1
>>> apart((x + 1)/(x - 1), x)
2
1 + -----
x - 1
将分数组合在一起,使用函数 together(expr, x) :
>>> from sympy import together
>>> together(1/x + 1/y + 1/z)
x*y + x*z + y*z
---------------
x*y*z
>>> together(apart((x + 1)/(x - 1), x), x)
x + 1
-----
x - 1
>>> together(apart(1/( (x + 2)*(x + 1) ), x), x)
1
---------------
(x + 1)*(x + 2)
在SymPy中使用极限很简单,使用函数 limit(function, variable, point) ,所以为了计算当 x -> 0 时f(x)的极限,只需要输入 limit(f, x, 0) :
>>> from sympy import limit, Symbol, sin, oo
>>> x = Symbol("x")
>>> limit(sin(x)/x, x, 0)
1
你还可以计算在无穷远处的极限:
>>> limit(x, x, oo)
oo
>>> limit(1/x, x, oo)
0
>>> limit(x**x, x, 0)
1
对于一些关于极限的非平凡的例子,你可以查阅测试文件 test_demidovich.py
你可以对任何SymPy的表达式进行求导,使用函数 diff(func, var) ,例如:
>>> from sympy import diff, Symbol, sin, tan
>>> x = Symbol('x')
>>> diff(sin(x), x)
cos(x)
>>> diff(sin(2*x), x)
2*cos(2*x)
>>> diff(tan(x), x)
2
tan (x) + 1
你可以检查,但是这是正确的:
>>> from sympy import limit
>>> from sympy.abc import delta
>>> limit((tan(x + delta) - tan(x))/delta, delta, 0)
2
tan (x) + 1
计算更高阶的导数,使用函数 diff(func, var, n) :
>>> diff(sin(2*x), x, 1)
2*cos(2*x)
>>> diff(sin(2*x), x, 2)
-4*sin(2*x)
>>> diff(sin(2*x), x, 3)
-8*cos(2*x)
使用 .series(var, point, order) :
>>> from sympy import Symbol, cos
>>> x = Symbol('x')
>>> cos(x).series(x, 0, 10)
2 4 6 8
x x x x / 10\
1 - -- + -- - --- + ----- + O\x /
2 24 720 40320
>>> (1/cos(x)).series(x, 0, 10)
2 4 6 8
x 5*x 61*x 277*x / 10\
1 + -- + ---- + ----- + ------ + O\x /
2 24 720 8064
另一个简单的例子:
>>> from sympy import Integral, pprint
>>> y = Symbol("y")
>>> e = 1/(x + y)
>>> s = e.series(x, 0, 5)
>>> print(s)
1/y - x/y**2 + x**2/y**3 - x**3/y**4 + x**4/y**5 + O(x**5)
>>> pprint(s)
2 3 4
1 x x x x / 5\
- - -- + -- - -- + -- + O\x /
y 2 3 4 5
y y y y
计算求和变量在给定的范围内函数f的和。
函数summation(f, (i, a, b))计算当i在[a, b]内时f的和是多少,即:
b
____
\ `
summation(f, (i, a, b)) = ) f
/___,
i = a
如果不能求和,会输出相应的求和公式。通过引入变量的额外限制,可以复合求和:
>>> from sympy import summation, oo, symbols, log
>>> i, n, m = symbols('i n m', integer=True)
>>> summation(2*i - 1, (i, 1, n))
2
n
>>> summation(1/2**i, (i, 0, oo))
2
>>> summation(1/log(n)**n, (n, 2, oo))
oo
___
\ `
\ -n
/ log (n)
/__,
n = 2
>>> summation(i, (i, 0, n), (n, 0, m))
3 2
m m m
-- + -- + -
6 2 3
>>> from sympy.abc import x
>>> from sympy import factorial
>>> summation(x**n/factorial(n), (n, 0, oo))
x
e
SymPy支持初等超越函数和特使函数的不定积分和定积分,通过使用函数 integrate() ,这个函数用到了强大的扩展Risch-Norman算法、一些启发式算法和模式匹配:
>>> from sympy import integrate, erf, exp, sin, log, oo, pi, sinh, symbols
>>> x, y = symbols('x,y')
求初等函数的积分:
>>> integrate(6*x**5, x)
6
x
>>> integrate(sin(x), x)
-cos(x)
>>> integrate(log(x), x)
x*log(x) - x
>>> integrate(2*x + sinh(x), x)
2
x + cosh(x)
求特殊函数的积分也很容易:
>>> integrate(exp(-x**2)*erf(x), x)
____ 2
\/ pi *erf (x)
--------------
4
计算定积分:
>>> integrate(x**3, (x, -1, 1))
0
>>> integrate(sin(x), (x, 0, pi/2))
1
>>> integrate(cos(x), (x, -pi/2, pi/2))
2
SymPy还支持广义积分:
>>> integrate(exp(-x), (x, 0, oo))
1
>>> integrate(log(x), (x, 0, 1))
-1
除了虚数单位I,其它符号的创建都可以带有属性(例如实数、正数、复数等等),而这些属性会影响操作结果:
>>> from sympy import Symbol, exp, I
>>> x = Symbol("x") # a plain x with no attributes
>>> exp(I*x).expand()
I*x
e
>>> exp(I*x).expand(complex=True)
-im(x) -im(x)
I*e *sin(re(x)) + e *cos(re(x))
>>> x = Symbol("x", real=True)
>>> exp(I*x).expand(complex=True)
I*sin(x) + cos(x)
三角函数
>>> from sympy import asin, asinh, cos, sin, sinh, symbols, I
>>> x, y = symbols('x,y')
>>> sin(x + y).expand(trig=True)
sin(x)*cos(y) + sin(y)*cos(x)
>>> cos(x + y).expand(trig=True)
-sin(x)*sin(y) + cos(x)*cos(y)
>>> sin(I*x)
I*sinh(x)
>>> sinh(I*x)
I*sin(x)
>>> asinh(I)
I*pi
----
2
>>> asinh(I*x)
I*asin(x)
>>> sin(x).series(x, 0, 10)
3 5 7 9
x x x x / 10\
x - -- + --- - ---- + ------ + O\x /
6 120 5040 362880
>>> sinh(x).series(x, 0, 10)
3 5 7 9
x x x x / 10\
x + -- + --- + ---- + ------ + O\x /
6 120 5040 362880
>>> asin(x).series(x, 0, 10)
3 5 7 9
x 3*x 5*x 35*x / 10\
x + -- + ---- + ---- + ----- + O\x /
6 40 112 1152
>>> asinh(x).series(x, 0, 10)
3 5 7 9
x 3*x 5*x 35*x / 10\
x - -- + ---- - ---- + ----- + O\x /
6 40 112 1152
球谐函数
>>> from sympy import Ylm
>>> from sympy.abc import theta, phi
>>> Ylm(1, 0, theta, phi)
___
\/ 3 *cos(theta)
----------------
____
2*\/ pi
>>> Ylm(1, 1, theta, phi)
___ I*phi
-\/ 6 *e *sin(theta)
------------------------
____
4*\/ pi
>>> Ylm(2, 1, theta, phi)
____ I*phi
-\/ 30 *e *sin(theta)*cos(theta)
------------------------------------
____
4*\/ pi
阶乘和gamma函数
>>> from sympy import factorial, gamma, Symbol
>>> x = Symbol("x")
>>> n = Symbol("n", integer=True)
>>> factorial(x)
x!
>>> factorial(n)
n!
>>> gamma(x + 1).series(x, 0, 3) # i.e. factorial(x)
/ 2 2\
2 |EulerGamma pi | / 3\
1 - EulerGamma*x + x *|----------- + ---| + O\x /
\ 2 12/
zeta函数
>>> from sympy import zeta
>>> zeta(4, x)
zeta(4, x)
>>> zeta(4, 1)
4
pi
---
90
>>> zeta(4, 2)
4
pi
-1 + ---
90
>>> zeta(4, 3)
4
17 pi
- -- + ---
16 90
多项式
>>> from sympy import assoc_legendre, chebyshevt, legendre, hermite
>>> chebyshevt(2, x)
2
2*x - 1
>>> chebyshevt(4, x)
4 2
8*x - 8*x + 1
>>> legendre(2, x)
2
3*x 1
---- - -
2 2
>>> legendre(8, x)
8 6 4 2
6435*x 3003*x 3465*x 315*x 35
------- - ------- + ------- - ------ + ---
128 32 64 32 128
>>> assoc_legendre(2, 1, x)
__________
/ 2
-3*x*\/ - x + 1
>>> assoc_legendre(2, 2, x)
2
- 3*x + 3
>>> hermite(3, x)
3
8*x - 12*x
在 isympy :
>>> from sympy import Function, Symbol, dsolve
>>> f = Function('f')
>>> x = Symbol('x')
>>> f(x).diff(x, x) + f(x)
2
d
f(x) + ---(f(x))
2
dx
>>> dsolve(f(x).diff(x, x) + f(x), f(x))
f(x) = C1*sin(x) + C2*cos(x)
在 isympy :
>>> from sympy import solve, symbols
>>> x, y = symbols('x,y')
>>> solve(x**4 - 1, x)
[-1, 1, -I, I]
>>> solve([x + 5*y - 2, -3*x + 6*y - 15], [x, y])
{x: -3, y: 1}
新建的矩阵是矩阵类的实例:
>>> from sympy import Matrix, Symbol
>>> Matrix([[1, 0], [0, 1]])
[1 0]
[ ]
[0 1]
矩阵还可以包含符号:
>>> x = Symbol('x')
>>> y = Symbol('y')
>>> A = Matrix([[1, x], [y, 1]])
>>> A
[1 x]
[ ]
[y 1]
>>> A**2
[x*y + 1 2*x ]
[ ]
[ 2*y x*y + 1]
想了解更多关于矩阵的内容,请看线性代数教程。
使用 .match() 模式和 Wild 类,可以对表达式执行模式匹配,返回结果就是所需要的各个匹配值,以字典的形式呈现,就像:
>>> from sympy import Symbol, Wild
>>> x = Symbol('x')
>>> p = Wild('p')
>>> (5*x**2).match(p*x**2)
{p: 5}
>>> q = Wild('q')
>>> (x**2).match(p*x**q)
{p: 1, q: 2}
如果匹配不成功,则返回 None :
>>> print (x + 1).match(p**x)
None
还可以在 Wild 类中添加排除参数,使得这些被排除的参数不会出现在结果当中:
>>> p = Wild('p', exclude=[1, x])
>>> print (x + 1).match(x + p) # 1 is excluded
None
>>> print (x + 1).match(p + 1) # x is excluded
None
>>> print (x + 1).match(x + 2 + p) # -1 is not excluded
{p_: -1}
打印表达式有很多方式。
标准打印
就是 str(expression) 的返回值,它看起来就像:
>>> from sympy import Integral
>>> from sympy.abc import x
>>> print x**2
x**2
>>> print 1/x
1/x
>>> print Integral(x**2, x)
Integral(x**2, x)
Pretty打印
使用 pprint 函数可以打印出好看的ASCII形式:
>>> from sympy import Integral, pprint
>>> from sympy.abc import x
>>> pprint(x**2)
2
x
>>> pprint(1/x)
1
-
x
>>> pprint(Integral(x**2, x))
/
|
| 2
| x dx
|
/
如果你安装了unicode字体,那么打印方式会默认使用 pprint 函数,可以使用 use_unicode 选项来替换这种打印方式:
>>> pprint(Integral(x**2, x), use_unicode=True)
⌠
⎮ 2
⎮ x dx
⌡
想获得更多unicode打印方式的例子,还可以查看维基百科 Pretty Printing 。
提示:通过这种方法可以在Python解释器中默认使用pretty打印:
$ python
Python 2.5.2 (r252:60911, Jun 25 2008, 17:58:32)
[GCC 4.3.1] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> from sympy import init_printing, var, Integral
>>> init_printing(use_unicode=False, wrap_line=False, no_global=True)
>>> var("x")
x
>>> x**3/3
3
x
--
3
>>> Integral(x**2, x) #doctest: +NORMALIZE_WHITESPACE
/
|
| 2
| x dx
|
/
Python打印
>>> from sympy.printing.python import python
>>> from sympy import Integral
>>> from sympy.abc import x
>>> print python(x**2)
x = Symbol('x')
e = x**2
>>> print python(1/x)
x = Symbol('x')
e = 1/x
>>> print python(Integral(x**2, x))
x = Symbol('x')
e = Integral(x**2, x)
LaTeX打印
>>> from sympy import Integral, latex
>>> from sympy.abc import x
>>> latex(x**2)
x^{2}
>>> latex(x**2, mode='inline')
$x^{2}$
>>> latex(x**2, mode='equation')
\begin{equation}x^{2}\end{equation}
>>> latex(x**2, mode='equation*')
\begin{equation*}x^{2}\end{equation*}
>>> latex(1/x)
\frac{1}{x}
>>> latex(Integral(x**2, x))
\int x^{2}\, dx
MathML
>>> from sympy.printing.mathml import mathml
>>> from sympy import Integral, latex
>>> from sympy.abc import x
>>> print mathml(x**2)
<apply><power/><ci>x</ci><cn>2</cn></apply>
>>> print mathml(1/x)
<apply><power/><ci>x</ci><cn>-1</cn></apply>
Pyglet
>>> from sympy import Integral, preview
>>> from sympy.abc import x
>>> preview(Integral(x**2, x))
如果安装了pyglet,一个呈现LaTeX表达式的窗口会被打开:
isympy 自动调用 pprint ,所以在默认情况下你看到的是pretty打印形式的输出。
值得注意的是,SymPy还提供一个打印模块 sympy.printing ,这个模块提供的打印方式有:
看完以上的教程之后,现在是时候进一步了解SymPy,请阅读 SymPy User’s Guide 和 SymPy Modules Reference 。
还可以浏览我们的公共网页 wiki.sympy.org ,上面有很多有用的例子、教程和cookbooks供大家参考,这些资料都是我们和用户编辑的,欢迎您的参与。
此教程还有其它的语言版本: