Python中用于符号型数学计算的sympy模块

Python sympy模块

SymPy一个用于符号型数学计算（symbolic mathematics）的Python库。它旨在成为一个功能齐全的计算机代数系统（Computer Algebra System，CAS），同时保持代码简洁、易于理解和扩展。SymPy完全是用Python写的，并不需要外部的库。

基础语法

符号表示

SymPy库中使用Symbol函数定义符号变量

from sympy import *
x=Symbol('x')
y=Symbol('y')
# ------------------
x, y=Symbol('x y') # 第二个用空格隔开

方程表示

使用代码表示数学符号与手写体的数学运算符号存在一定的差异，对于长的表达式，如果不确定运算符的优先级，可以加入小括号提升其优先级。下面列举常用的运算符：

• 加号 +

• 减号 -

• 除号 /

• 乘号 *

• 指数 **

• 对数 log()

• e的指数次幂 exp()

符号的输出设置

在sympy里进行符号运算之前，必须先定义sympy的符号，这样sympy才能识别该符号。

.init_printing(use_latex=True)开启时，运行时输出的是LaTeX的格式。

使用：latex()函数，同样返回LaTeX的格式。

import sympy
x, y, z = sympy.Symbol('x y z')  # 符号化变量
sympy.init_printing(use_latex=True)  # 输出设置
print("x:", type(x))
print("y:", type(y))
print(x ** 2 + y + z)
print(sympy.latex(x ** 2 + y + z))
# output
# x: <class 'sympy.core.symbol.Symbol'>
# y: <class 'sympy.core.symbol.Symbol'>
# x**2 + y + z
# x^{2} + y + z

替换符号

sub是Substitution的简称，也就是替换。语法是：expr.sub(old, new)其有两个作用：

1. 数值替换，用数值替换符号，进行带入计算。

2. 符号替换，用一些符号替换符号。

import sympy

x, y, z = sympy.symbols('x y z')
expr = x ** 2 + 1

# 数值替换
result = expr.subs(x, 2)
print("原式：", expr)
print("数值计算：", result)

# 符号替换
new_expr = expr.subs(x, y + z)
print("符号替换：", new_expr)

# output
# 原式： x**2 + 1
# 数值计算： 5
# 符号替换： (y + z)**2 + 1

PS：

• subs()函数不改变原表达式，并且返回一个修改的表达式。

• 当需要替换多个表达式时，可以在subs()里使用列表。如：subs([(x,2), (y, 3), (z, 4)])，表示：将x替换成2，y替换成3，z替换成4。

将字符串变为sympy的表达式

import sympy
string = "x**2+2*y + z/2"
expr = sympy.sympify(string)  # 转化
print("类型：", type(expr))
print("表达式：", expr)
# output
# 类型： <class 'sympy.core.add.Add'>
# 表达式： x**2 + 2*y + z/2

PS：不要混淆了sympify()函数与 simplify()函数，前者是转化，后者是简化。

数值计算

相当于python自带的eval()函数，只是进行的是float浮点数运算。

（1）对于数值表达式的运算：

import sympy
expr = sympy.sqrt(8)
result = expr.evalf()  # 进行计算
print(result)
# output
# 2.82842712474619

（2）对于符号表达式的运算：

对于表达式常用的是:.evalf(subs={x: 2.4})

import sympy
x = sympy.Symbol('x')
expr = x**2+3
result = expr.evalf(subs={x: 2})
print(result)
# output
# 7.00000000000000

自定义表达式

该函数有点类似于lambda()，用于自己构造一个函数表达。

import sympy
import numpy as np
x = sympy.Symbol('x')
a = np.arange(10)
expr = x ** 2
f = sympy.lambdify(x, expr, "numpy")  # 构造自己的函数
print(f(a))
# output
# [ 0  1  4  9 16 25 36 49 64 81]

解方程

解方程的功能主要是使用Sympy中solve函数实现。Solve函数的第一个参数是要解的方程，要求右端等于0，第二个参数是未知数。

import sympy
x = sympy.symbols("x")  # 申明未知数"x"
a = sympy.solve([x + (1 / 5) * x - 240], [x])  # 写入需要解的方程体
print(a) 

在写入方程的时候，将等号右边的数移到了等号左边，然后将等号丢弃，最后变成了一个式子。（注意移动过程中的变号）

解方程组：

import sympy
x, y = sympy.symbols("x y")
a = sympy.solve([3 * x - 2 * y - 3, x + 2 * y - 5], [x, y])
print(a)

解分式方程：

import sympy
x, y = sympy.symbols("x y")
a = sympy.solve([((x + 1) / x + 1 / (x - 2)) - 1], [x])
print(a)

巧算代数式，比如下面这道题：

import sympy
x, y = sympy.symbols("x y")
a = sympy.solve([x + y - 0.2, x + 3 * y - 1], [x, y])
x, y = a[x], a[y]
re = x ** 2 + 4 * x * y + 4 * y ** 2
print(re)

求极限

需要用到limit()函数求极限。

import sympy
n = sympy.symbols("n")
s = ((n + 3) / (n + 2)) ** n
a = sympy.limit(s, n, "oo") 
print(a)

无穷的表示方法是两个小写的字母o，即"oo"或"-oo"。

求不定积分

用法其实和上面都一样。

import sympy
x = sympy.symbols("x")
a = sympy.integrate("x**2", x)
print(a)

integrate是不定积分运算函数。这个函数有两个参数：

• 前面一个参数是需要计算不定积分的函数

• 后面的参数是以谁为参数计算不定积分。

求定积分

import sympy
t = sympy.Symbol("t")
x = sympy.Symbol("x")
m = sympy.integrate("sin(t) / (pi - t)", (t, 0, x))
n = sympy.integrate(m, (x, 0, "pi"))
print(n)

求解微分

diff函数：diff(f(x), x, k)。k表示求k阶导的意思。

import sympy
f = sympy.Function("f")
x = sympy.symbols("x")
a = sympy.diff(x**3, x, 2)
print(a)

解微分方程

dsolve函数是用来解决微分方程的函数。函数的一个用法为：dsolve(eq, f(x))。

• 第一个参数为微分方程（要先将等式移项为右端为0的形式)

• 第二个参数为要解的函数(在微分方程中)

import sympy
f = sympy.Function("f")
x = sympy.symbols("x")
a = 2 * x - sympy.diff(f(x), x)
b = sympy.dsolve(a, f(x))
print(b)  # Eq(f(x), C1 + x**2)

最后的答案Eq(f(x), C1*exp(x**2))表示：f(x) = C1*exp(x**2)

矩阵化简

import sympy
x1, x2, x3 = sympy.symbols('x1 x2 x3')
a11, a12, a13, a22, a23, a33 = sympy.symbols('a11 a12 a13 a22 a23 a33')
m = sympy.Matrix([[x1, x2, x3]])  # 注意有两个括号
n = sympy.Matrix([[a11, a12, a13], [a12, a22, a23], [a13, a23, a33]])
v = sympy.Matrix([[x1], [x2], [x3]])
print(m * n * v)  # Matrix([[x1*(a11*x1 + a12*x2 + a13*x3) + x2*(a12*x1 + a22*x2 + a23*x3) + x3*(a13*x1 + a23*x2 + a33*x3)]])

如果上式中x1，x2，x3均等于1，则可这样代入：

import sympy
x1, x2, x3 = sympy.symbols('x1 x2 x3')
a11, a12, a13, a22, a23, a33 = sympy.symbols('a11 a12 a13 a22 a23 a33')
m = sympy.Matrix([[x1, x2, x3]])  # 注意有两个括号
n = sympy.Matrix([[a11, a12, a13], [a12, a22, a23], [a13, a23, a33]])
v = sympy.Matrix([[x1], [x2], [x3]])
f = m * n * v
print(f[0].subs({x1: 1, x2: 1, x3: 1}))  # a11 + 2*a12 + 2*a13 + a22 + 2*a23 + a33