SymPy实现高等数学符号解

一些前置知识:

x=symbols('x')声明变量\(x\)

y.subs(x,x0)\(y\)\(x\)的表达式时,带入\(x=x_0\)的值

y.n()常量表达式求浮点值

求极限

limit(e, z, z0, dir="+")

\(e\)\(z\)的表达式,求\(e(z)\)\(z_0\)处的极限,\(\text{dir}\)可以指定是\(z_0^+\)还是\(z_0^-\)

以下计算了\(\lim_{x\rightarrow 0}\frac{\sin x}x\)\(\lim_{x\rightarrow \infty} (1+\frac 1 x)^x\)\(\lim_{x\rightarrow 0}\sin(\frac 1 x)\)三个极限

1
2
3
4
5
from sympy import *
x = symbols('x')
print(limit(sin(x)/x, x, 0))
print(limit(pow(1+1/x, x), x, oo))
print(limit(sin(1/x), x, 0))

输出:

1
2
3
1
E
AccumBounds(-1, 1)

求导数

diff(expr,x1,n1,x2,n2,...)

表达式expr依次对xini阶偏导,ni不写默认为1

\(z=\sin x+x^2e^y\),求\(\frac{\partial^2 z}{\partial x^2},\frac{\partial z}{\partial y},\frac{\partial^2 z}{\partial x\partial y}\)

1
2
3
4
5
6
from sympy import *
x, y = symbols("x y")
z = sin(x)+x**2*exp(y)
print(diff(z, x, 2))
print(diff(z, y))
print(diff(z, x, y))

输出:

1
2
3
2*exp(y) - sin(x)
x**2*exp(y)
2*x*exp(y)

级数求和

summation(expr,(i,a,b))

expr是关于i的函数,求\(\sum_{i=a}^b expr(i)\)

示例:\(\sum_{k=1}^n k^2\)\(\sum_{k=1}^\infty \frac 1 {k^2}\)

1
2
3
4
from sympy import *
k, n = symbols('k n')
print(factor(summation(k**2, (k, 1, n)))) # 计算结果因式分解
print(summation(1/k**2, (k, 1, oo)))

输出:

1
2
n*(n + 1)*(2*n + 1)/6
pi**2/6

泰勒展开

series(expr, x=None, x0=0, n=6, dir="+")

\(expr\)\(x\)的表达式,在\(x_0\)处的\(n\)阶展开,\(\text{dir}\)可以指定是\(x_0^+\)还是\(x_0^-\)

下面演示了\(\sin x\)在0点的3,5,7阶泰勒展开式

1
2
3
4
5
6
7
8
9
10
11
12
from pylab import rc
from sympy import *
from sympy.plotting import *
rc('font', size=16)
x = symbols('x')
y = sin(x)
yi = list(range(8))
for k in range(3, 8, 2):
print(series(y, x, 0, k))
yi[k] = series(y, x, 0, k).removeO() # 带有余项的表达式无法作图
plot(y, yi[3], yi[5], yi[7], (x, 0, 2),
xlabel='$x$', ylabel='$y$').save('series.png')

输出:

1
2
3
x + O(x**3)
x - x**3/6 + O(x**5)
x - x**3/6 + x**5/120 + O(x**7)

积分

不定积分:

integrate(f,x)

相当于\(\int f(x)\mathrm d x\)

定积分:

integrate(f,(x,a,b))

相当于\(\int_a^b f(x)\mathrm d x\)

以下代码求解\(\int_0^\pi \sin(2x)\mathrm d x,\int_0^{+\infty} \frac {\sin(x)} x\mathrm d x,\int \sin(x)e^x\mathrm d x\)

1
2
3
4
5
from sympy import *
x = symbols('x')
print(integrate(sin(2*x), (x, 0, pi)))
print(integrate(sin(x)/x, (x, 0, oo)))
print(integrate(sin(x)*exp(x), x))

输出:

1
2
3
0
pi/2
exp(x)*sin(x)/2 - exp(x)*cos(x)/2

解代数方程(组)

solve(expr,x)

求解方程\(expr=0\),其中\(x\)是自变量

返回一个list存放所有的根

solve([...,expr_i,...],[...,x_i,...])

求解方程组,其中\(expr_i\)对应\(x_i\)

返回一个list存放所有根,每个根用一个tuple表示,对应此根所有各个变量的值

roots(expr,x)

返回一个dict,表示每个根以及它的重数

示例: \[ x^3-1=0 \\ (x-2)^2(x-1)^3=0 \\ \begin{equation} \begin{cases} x^2+y^2=1 \\ x-y=0 \end{cases} \end{equation} \]

1
2
3
4
5
6
from sympy import *
x, y = symbols("x y")
print(solve(x**3-1, x))
print(solve((x-2)**2*(x-1)**3, x))
print(roots((x-2)**2*(x-1)**3, x))
print(solve([x**2+y**2-1, x-y], [x, y]))

输出:

1
2
3
4
[1, -1/2 - sqrt(3)*I/2, -1/2 + sqrt(3)*I/2]
[1, 2]
{2: 2, 1: 3}
[(-sqrt(2)/2, -sqrt(2)/2), (sqrt(2)/2, sqrt(2)/2)]

微分方程(组)的通解

微分方程中的\(y\)是关于\(x\)的函数,而且是未知函数,而非未知量,因此要将\(y\)声明为Function

y=Function('y')

y=symbols('y',cls=Function)

两种写法都是可以的

dsolve(eq,f(x),ics=),其中eq是关于f(x)的表达式,求解\(eq=0\)

ics可以以一个dict指定f(x)的初值

示例: \[ y''-5y'+6y=0 \\ y''-5y'+6y=xe^{2x} \\ y''-5y'+6y=0,y(0)=1,y'(0)=0 \\ y''-5y'+6y=xe^{2x},y(0)=1,y(2)=0 \]

1
2
3
4
5
6
7
8
9
from sympy import *
x = symbols('x')
y = Function('y')
eq1 = diff(y(x), x, 2)-5*diff(y(x), x)+6*y(x)
eq2 = diff(y(x), x, 2)-5*diff(y(x), x)+6*y(x)-x*exp(2*x)
print(dsolve(eq1, y(x)))
print(dsolve(eq2, y(x)))
print(dsolve(eq1, y(x), ics={y(0): 1, diff(y(x), x).subs(x, 0): 0}))
print(dsolve(eq2, y(x), ics={y(0): 1, y(2): 0}))

输出:

1
2
3
4
Eq(y(x), (C1 + C2*exp(x))*exp(2*x))
Eq(y(x), (C1 + C2*exp(x) - x**2/2 - x)*exp(2*x))
Eq(y(x), (3 - 2*exp(x))*exp(2*x))
Eq(y(x), (-x**2/2 - x - 3*exp(x)/(1 - exp(2)) + (4 - exp(2))/(1 - exp(2)))*exp(2*x))