使用Python的SymPy库求解不定积分

如题。

参考

https://docs.sympy.org/latest/tutorials/intro-tutorial/intro.html
https://blog.csdn.net/Havingfunin/article/details/104741277
【sympy】用python的库 sympy 求积分
https://www.jianshu.com/p/3aea36f7a231
https://www.zhihu.com/question/414548947
用sympy求解积分的一个小例子 - 人生不会倒流的文章 - 知乎
SymPy 符号计算基本教程 - Longson的文章 - 知乎

计算机代数系统(Computer Algebra System),简称CAS,其标志是能够以字符串作为单位进行运算。科学计算可分为两类,一类是纯数值的计算,例如求函数的值,方程的数值解;另一类计算就是符号计算,这就像我们平时在数学的教学和研究中用笔和纸进行的数学运算一样,用各种字母 $x,y,z$ 代替具体数值,进行算式和方程的推导。

常见的编程语言如python、R、C++、Java等都支持数值计算,用户输入算式,程序输出结果。但是有些时候我们会碰到这样的任务:函数求导(不需要知道某个点的导函数,而是求解整个导函数的表达式)、不定积分(同样是求表达式)、多项式化简等。此时,数值计算就无能为力了。

SymPy是一个符号计算的Python库。它的目标是成为一个全功能的计算机代数系统,同时保持代码简洁、易于理解和扩展。它完全由Python写成,不依赖于外部库。

在下面这篇文章中,作者详细列出了sympy可以做的事情,本文不再赘述。

《SymPy 符号计算基本教程》 - Longson的文章 - 知乎

本文重点讲述,如何使用sympy计算不定积分:

Sympy是使用integrate(表达式,变量)来求不定积分的,具体来说:

  • integrate(f(x),x) 计算不定积分 $\int{f(x)dx}$
  • integrate(f(x),(x,a,b)) 计算定积分 $\int_a^b{f(x)dx}$

例如

求解定积分

$$
F(x)=\int_0^1x^{2} + e^{x} + 1\ \text{d}x
$$

1
2
3
from sympy import *
x = symbols('x')
print(integrate(x**2 + exp(x) + 1, (x, 0, 1))) # 定积分

output:

1
1/3 + E

上述程序如果求解的是不定积分,则代码改为

1
2
3
from sympy import *
x = symbols('x')
print(integrate(x**2 + exp(x) + 1, x)) # 定积分

output:

1
x**3/3 + x + exp(x)

即原函数为

$$
F(x)=\frac{x^{3}}{3} + x + e^{x}
$$

另一个例子:求解不定积分

$$
F(x)=\int\frac{x e^{x}}{\left(e^{x} + 1\right)^{2}} \text{d}x
$$

1
2
3
4
from sympy import *
x = symbols('x')
print(integrate(x*exp(x)/(exp(x)+1)**2, x))

output

1
x - x/(exp(x) + 1) - log(exp(x) + 1)

也就是说,原函数为

$$
F(x)=x - \frac{x}{e^{x} + 1} - \log{\left(e^{x} + 1 \right)}
$$

当然,在某些情况下,sympy可能没法一下子就求出原函数,此时求积分可能需要一些技巧:

例如求解定积分(来自知乎问题 《这个sympy 库为什么求不出这个函数的不定积分?》

$$
F(x)=\int_0^{\pi/2}\sqrt{1 - \cos{\left(2 x \right)}} \text{d}x
$$

1
2
3
4
5
6
7
import sympy as sym
from sympy import sin,cos,sec,exp,asin,sqrt,pi #必须将特殊的函数和常数导入程序才能运行
x = sym.symbols('x') #定义符号变量
y = sym.symbols('y')

d = sym.Integral(sqrt(1-cos(2*x)),(x,0,pi/2))
print(d)

输出依然为

1
Integral(sqrt(1 - cos(2*x)), (x, 0, pi/2))

也就是说,并没有求出具体的解。

但是我们可以用一些方法得到问题的解(参考文章 《用sympy求解积分的一个小例子 - 人生不会倒流的文章 - 知乎》

1
2
3
4
5
6
7
8
9
import sympy as sym
from sympy import sin,cos,sec,exp,asin,sqrt,pi #必须将特殊的函数和常数导入程序才能运行
x = sym.symbols('x') #定义符号变量
y = sym.symbols('y')

d = sym.Integral(sqrt(1-cos(2*x)),(x,0,pi/2))
print(d)
print(d.simplify())
print(d.simplify().simplify())

输出:

1
2
3
4
5
6
7
8
9
10
In [20]: d = sym.Integral(sqrt(1-cos(2*x)),(x,0,pi/2))

In [21]: print(d)
Integral(sqrt(1 - cos(2*x)), (x, 0, pi/2))

In [22]: print(d.simplify())
Integral(sqrt(2)*sqrt(sin(x)**2), (x, 0, pi/2))

In [23]: print(d.simplify().simplify())
sqrt(2)

integrate()的结果,第一次使用simplify()进行化简,得到化简后的表达式,此时的表达式已经能够进行求值了。

$$
\int\limits_{0}^{\frac{\pi}{2}} \sqrt{2} \sqrt{\sin^{2}{\left(x \right)}}\ \text{d}x
$$

第二次使用simplify()进行化简,则直接得到所求定积分的结果 $\sqrt{2}$ 。这真的很神奇!