python log函數_用 Python 來研究數學 --- SymPy 符号工具包介紹

用 PYTHON 來研究數學 — SYMPY 符号工具包介紹

SymPy 的簡單介紹

SymPy

是一個

符号計算

的 Python 庫，完全由 Python 寫成，為許多數值分析，符号計算提供了重要的工具。SymPy 的第一個版本于 2007 年開源，并且經曆了十幾個版本的疊代，在 2019 年已經基于修正的 BSD 許可證開源了 1.4 版本。SymPy 的開源位址和官方網站分别是：

1. GitHub 連結：https://github.com/sympy/sympy
2. SymPy 官方網站：https://www.sympy.org/en/index.html

SymPy 的 logo

從

SymPy

的 1.4 版本文檔中，可以看出，SymPy 可以支援很多初等數學，高等數學，甚至研究所學生數學的符号計算。在初等數學和高等數學中，SymPy 可以支援的内容包括但不限于：

1. 基礎計算（Basic Operations）；
2. 公式簡化（Simplification）；
3. 微積分（Calculus）；
4. 解方程（Solver）；
5. 矩陣（Matrices）；
6. 幾何（geometry）；
7. 級數（Series）；

在更多的數學領域中，SymPy 可以支援的内容包括但不限于：

1. 範疇論（Category Theory）；
2. 微分幾何（Differential Geometry）；
3. 常微分方程（ODE）；
4. 偏微分方程（PDE）；
5. 傅立葉變換（Fourier Transform）；
6. 集合論（Set Theory）；
7. 邏輯計算（Logic Theory）。

SymPy 的教學目錄

SymPy 的工具庫介紹

SymPy 的基礎計算

在數學中，基礎的計算包括實數和複數的加減乘除，那麼就需要在程式中描述出實數與複數。著名的歐拉公式

。

正好用到了數學中最常見的五個實數。在 SymPy 裡面，

是用以下符号來表示的：其中 sympy.exp() 表示以

為底的函數。

sympy.exp(1), sympy.I, sympy.pi, sympy.oo
           

而想要計算歐拉公式的話，隻需要輸入下面的公式即可：

>>> sympy.exp(sympy.I * sympy.pi) + 1
0
           

如果需要看

的小數值，可以使用 evalf() 函數，其中 evalf() 函數裡面的值表示有效數字的位數。例如下面就是精确到 10 位有效數字。當然，也可以不輸入。

>>> sympy.E.evalf(10)
2.718281828
>>> sympy.E.evalf()
2.71828182845905
>>> sympy.pi.evalf(10)
3.141592654
>>> sympy.pi.evalf()
3.14159265358979
           

除此之外，如果需要檢視某個實數的有效數字，也是類似操作的：

>>> expr = sympy.sqrt(8)
>>> expr.evalf()
2.82842712474619
           

而對于實數的加減乘除，則可以如下操作：

>>> x, y= sympy.symbols("x y")
>>> x + y
x + y
>>> x - y
x - y
>>> x * y
x*y
>>> x / y
x/y
           

而對于複數的加減乘除，則是類似的操作，令兩個複數分别是

和

。

>>> x1, y1, x2, y2 = sympy.symbols("x1 y1 x2 y2")
>>> z1 = x1 + y1 * sympy.I
x1 + I*y1
>>>  z2 = x2 + y2 * sympy.I
x2 + I*y2
>>> z1 + z2
x1 + x2 + I*y1 + I*y2
>>> z1 - z2
x1 - x2 + I*y1 - I*y2
>>> z1 * z2
(x1 + I*y1)*(x2 + I*y2)
>>> z1 / z2
(x1 + I*y1)/(x2 + I*y2)
           

對于多項式而言，有的時候我們希望将其展開，有的時候則需要将其合并，最終将其簡化成最簡單的形式。

>>> sympy.expand((x+1)**2)
x**2 + 2*x + 1
>>> sympy.expand((x+1)**5)
x**5 + 5*x**4 + 10*x**3 + 10*x**2 + 5*x + 1
>>> sympy.factor(x**3+1)
(x + 1)*(x**2 - x + 1)
>>> sympy.factor(x**2+3*x+2)
(x + 1)*(x + 2)
>>> sympy.simplify(x**2 + x + 1 - x)
x**2 + 1
>>> sympy.simplify(sympy.sin(x)**2 + sympy.cos(x)**2)
1
           

在多變量的場景下，SymPy 也可以對其中的某個變量合并同類項，同時還可以計算某個變量的某個次數所對應的系數是多少，例如：

>>> expr = x*y + x - 3 + 2*x**2 - x**2 + x**3 + y**2 + x**2*y**2
>>> sympy.collect(expr,x)
x**3 + x**2*(y**2 + 1) + x*(y + 1) + y**2 - 3
>>> sympy.collect(expr,y)
x**3 + x**2 + x*y + x + y**2*(x**2 + 1) - 3
>>> expr.coeff(x, 2)
y**2 + 1
>>> expr.coeff(y, 1)
x
           

有理函數形如

，其中

和

都是多項式。一般情況下，我們希望對有理函數進行簡化，合并或者分解的數學計算。

在需要合并的情形下，如果想把有理函數處理成标準格式

并且去除公因子，那麼可以使用 cancel 函數。另一個類似的就是 together 函數，但是不同之處在于 cancel 會消除公因子，together 不會消除公因子。例如：

，

>>> expr = (x**2 + 3*x + 2)/(x**2 + x)
>>> sympy.cancel(expr)
(x + 2)/x
>>> sympy.together(expr)
(x**2 + 3*x + 2)/(x*(x + 1))
           

除了合并和消除公因子之外，有的時候還希望對分子和分母進行因式分解，例如：

expr = (x**2 + 3*x + 2)/(x**2 + x)
>>> sympy.factor(expr)
(x + 2)/x
>>> expr = (x**3 + 3*x**2 + 2*x)/(x**5+x)
>>> sympy.factor(expr)
(x + 1)*(x + 2)/(x**4 + 1)
>>> expr = x**2 + (2*x+1)/(x**3+1)
>>> sympy.factor(expr)
(x**5 + x**2 + 2*x + 1)/((x + 1)*(x**2 - x + 1))
           

合并的反面就是部分分式展開（Partial Fraction Decomposition），它是把有理函數分解成多個次數較低的有理函數和的形式。這裡需要用 apart 函數：

>>> expr = (x**4 + 3*x**2 + 2*x)/(x**2+x)
>>> sympy.apart(expr)
x**2 - x + 4 - 2/(x + 1)
>>> expr = (x**5 + 1)/(x**3+1)
>>> sympy.apart(expr)
x**2 - (x - 1)/(x**2 - x + 1)
           

在 SymPy 裡面，同樣支援各種各樣的三角函數，例如：三角函數的簡化函數 trigsimp，三角函數的展開 expand_trig，

>>> expr = sympy.sin(x)**2 + sympy.cos(x)**2
>>> sympy.trigsimp(expr)
1
>>> sympy.expand_trig(sympy.sin(x+y))
sin(x)*cos(y) + sin(y)*cos(x)
>>> sympy.expand_trig(sympy.cos(x+y))
-sin(x)*sin(y) + cos(x)*cos(y)
>>> sympy.trigsimp(sympy.sin(x)*sympy.cos(y) + sympy.sin(y)*sympy.cos(x))
sin(x + y)
>>> sympy.trigsimp(-sympy.sin(x)*sympy.sin(y) + sympy.cos(x)*sympy.cos(y))
cos(x + y)
           

同樣的，在乘幂上面，同樣有簡化函數 powsimp，效果與之前提到的 simplify 一樣。除此之外，還可以根據底數來做合并，即分别使用 expand_power_exp 函數與 expand_power_base 函數。

>>> sympy.powsimp(x**z*y**z*x**z)
x**(2*z)*y**z
>>> sympy.simplify(x**z*y**z*x**z)
x**(2*z)*y**z
>>> sympy.expand_power_exp(x**(y + z))
x**y*x**z
>>> sympy.expand_power_base(x**(y + z))
x**(y + z)
           

作為指數的反函數對數，sympy 也是有着類似的展開合并函數，expand_log，logcombine 承擔了這樣的角色。

，

。

>>> sympy.expand_log(sympy.log(x*y), force=True)
log(x) + log(y)
>>> sympy.expand_log(sympy.log(x/y), force=True)
log(x) - log(y)
           

SymPy 的微積分工具

下面，我們會從一個最基本的函數

出發，來介紹 SymPy 的各種函數使用方法。如果想進行變量替換，例如把

變成

，那麼可以使用 substitution 方法。除此之外，有的時候也希望能夠得到函數

在某個點的取值，例如

，那麼可以把參數換成 1 即可得到函數的取值。例如，

>>> import sympy
>>> x = sympy.Symbol("x")
>>> f = 1 / x
1/x
>>> y = sympy.Symbol("y")
>>> f = f.subs(x,y)
1/y
>>> f = f.subs(y,1)
1
           

在微積分裡面，最常見的概念就是極限，SymPy 裡面的極限函數就是 limit。使用方法如下：

>>> f = 1/x
>>> sympy.limit(f,x,0)
oo
>>> sympy.limit(f,x,2)
1/2
>>> sympy.limit(f,x,sympy.oo)
0
>>> g = x * sympy.log(x)
>>> sympy.limit(g,x,0)
0
           

對于函數

而言，它的導數計算函數是 diff，

階導數也可以用這個函數算。

>>> f = 1/x
>>> sympy.diff(f,x)
-1/x**2
>>> sympy.diff(f,x,2)
2/x**3
>>> sympy.diff(f,x,3)
-6/x**4
>>> sympy.diff(f,x,4)
24/x**5
           

提到

階導數，就必須要提一下 Taylor Series 了。對于常見函數的 Taylor Series，SymPy 也是有非常簡便的方法，那就是 series 函數。其參數包括 expr, x, x0, n, dir，分别對應着表達式，函數的自變量，Taylor Series 的中心點，n 表示階數，dir 表示方向，包括"+-"，"-"，"+"，分别表示

,

,

。

sympy.series.series.series(expr,?x=None,?x0=0,?n=6,?dir='+')
>>> g = sympy.cos(x)
>>> sympy.series(g, x) 1 - x**2/2 + x**4/24 + O(x**6)
>>> sympy.series(g, x, x0=1, n=10) cos(1) - (x - 1)*sin(1) - (x - 1)**2*cos(1)/2 + (x - 1)**3*sin(1)/6 + (x - 1)**4*cos(1)/24 - (x - 1)**5*sin(1)/120 - (x - 1)**6*cos(1)/720 + (x - 1)**7*sin(1)/5040 + (x - 1)**8*cos(1)/40320 - (x - 1)**9*sin(1)/362880 + O((x - 1)**10, (x, 1))
           

積分的計算函數是 integrate，包括定積分與不定積分：

，

。

>>> f = 1/x
>>> sympy.integrate(f,x)
log(x)
>>> sympy.integrate(f, (x,1,2))
log(2)
           

對于廣義積分而言，就需要用到

這個概念了，但是在 SymPy 裡面的寫法還是一樣的。

，

，

。

>>> g = sympy.exp(-x**2)
>>> sympy.integrate(g, (x,-sympy.oo,0))
sqrt(pi)/2
>>> g = sympy.exp(-x)
>>> sympy.integrate(g, (x, 0, sympy.oo))
1
>>> h = sympy.exp(-x**2 - y**2)
>>> sympy.integrate(h, (x,-sympy.oo, sympy.oo), (y, -sympy.oo, sympy.oo))
pi
           

SymPy 的方程工具

在初等數學中，通常都存在一進制一次方程，一進制二次方程等，并且在不同的域上有着不同的解。SymPy 裡面的相應函數就是 solveset，根據定義域的不同，可以獲得完全不同的解。

，

，

，

，

。

>>> sympy.solveset(sympy.Eq(x**3,1), x, domain=sympy.S.Reals)
{1}
>>> sympy.solveset(sympy.Eq(x**3,1), x, domain=sympy.S.Complexes)
{1, -1/2 - sqrt(3)*I/2, -1/2 + sqrt(3)*I/2}
>>> sympy.solveset(sympy.Eq(x**3 - 1,0), x, domain=sympy.S.Reals)
{1}
>>> sympy.solveset(sympy.Eq(x**3 - 1,0), x, domain=sympy.S.Complexes)
{1, -1/2 - sqrt(3)*I/2, -1/2 + sqrt(3)*I/2}
>>> sympy.solveset(sympy.exp(x),x)
EmptySet()
>>> sympy.solveset(sympy.exp(x)-1,x,domain=sympy.S.Reals)
{0}
>>> sympy.solveset(sympy.exp(x)-1,x,domain=sympy.S.Complexes)
ImageSet(Lambda(_n, 2*_n*I*pi), Integers)
           

在這裡，Lambda 表示的是數學公式，第一個是自變量，第二個是函數，最後是自變量的定義域。

線上性代數中，最常見的還是多元一次方程組，那麼解法是一樣的：

,

>>> sympy.solve([x+y-10, x-y-2], [x,y])
{x: 6, y: 4}
           

對于三角函數，也是類似的寫法：

，

>>> sympy.solve([sympy.sin(x-y), sympy.cos(x+y)], [x,y])
[(-pi/4, 3*pi/4), (pi/4, pi/4), (3*pi/4, 3*pi/4), (5*pi/4, pi/4)]
           

SymPy 的矩陣工具

在矩陣論中，最常見的就是機關矩陣了，而機關矩陣隻與一個參數有關，那就是矩陣的大小。下面就是 3*3，3*2，2*3 大小的矩陣。

>>> sympy.eye(3)
Matrix([
[1, 0, 0],
[0, 1, 0],
[0, 0, 1]])
>>> sympy.eye(3,2)
Matrix([
[1, 0],
[0, 1],
[0, 0]])
>>> sympy.eye(2,3)
Matrix([
[1, 0, 0],
[0, 1, 0]])
           

另外還有全零和全一矩陣，意思就是矩陣中的所有值全部是零和一。

>>> sympy.ones(2,3)
Matrix([
[1, 1, 1],
[1, 1, 1]])
>>> sympy.zeros(3,2)
Matrix([
[0, 0],
[0, 0],
[0, 0]])
           

而對角矩陣也可以使用 diag 輕松獲得：

>>> sympy.diag(1,1,2)
Matrix([
[1, 0, 0],
[0, 1, 0],
[0, 0, 2]])
           

而矩陣的加法，減法，乘法，逆運算，轉置，行列式，SymPy 都是可以支援的：

，

，

>>> A = sympy.Matrix([[1,1],[0,2]])
>>> B = sympy.Matrix([[1,0],[1,1]])
>>> A
Matrix([
[1, 1],
[0, 2]])
>>> B
Matrix([
[1, 0],
[1, 1]])
>>> A + B
Matrix([
[2, 1],
[1, 3]])
>>> A - B
Matrix([
[ 0, 1],
[-1, 1]])
>>> A * B
Matrix([
[2, 1],
[2, 2]])
>>> A * B**-1
Matrix([
[ 0, 1],
[-2, 2]])
>>> B**-1
Matrix([
[ 1, 0],
[-1, 1]])
>>> A.T
Matrix([
[1, 0],
[1, 2]])
>>> A.det()
2
           

在某些情況下，需要對矩陣進行加上一行或者加上一列的操作，在這裡就可以使用 row_insert 或者 col_insert 函數：第一個參數表示插入的位置，第二個參數就是相應的行向量或者列向量。而在删除上就很簡單了，直接使用 row_del 或者 col_del 即可。

>>> A
Matrix([
[1, 1],
[0, 2]])
>>> A.row_insert(1, sympy.Matrix([[10,10]]))
Matrix([
[ 1, 1],
[10, 10],
[ 0, 2]])
>>> A.col_insert(0, sympy.Matrix([3,3]))
Matrix([
[3, 1, 1],
[3, 0, 2]])
>>> A.row_del(0)
>>> A
Matrix([[0, 2]])
>>> A.col_del(1)
>>> A
Matrix([[0]])
           

在對角化方面，同樣可以使用 eigenvals()，eigenvecs(), diagonalize() 函數：

>>> A
Matrix([
[1, 1],
[0, 2]])
>>> A.eigenvals()
{2: 1, 1: 1}
>>> A.eigenvects()
[(1, 1, [Matrix([
[1],
[0]])]), (2, 1, [Matrix([
[1],
[1]])])]
>>> A.diagonalize()
(Matrix([
[1, 1],
[0, 1]]), Matrix([
[1, 0],
[0, 2]]))
           

在 eigenvals() 傳回的結果中，第一個表示特征值，第二個表示該特征值的重數。在特征向量 eigenvecs() 中，第一個表示特征值，第二個表示特征值的重數，第三個表示特征向量。在對角化 diagonalize() 中，第一個矩陣表示

，第二個矩陣表示

，

。

在矩陣中，最常見的還是多元一次方程的解。如果要求

的解，可以有以下方案：

>>> A = sympy.Matrix([[1,1],[0,2]])
>>> A
Matrix([
[1, 1],
[0, 2]])
>>> b = sympy.Matrix([3,5])
>>> b
Matrix([
[3],
[5]])
>>> A**-1*b
Matrix([
[1/2],
[5/2]])
>>> sympy.linsolve((A,b))
{(1/2, 5/2)}
>>> sympy.linsolve((A,b),[x,y])
{(1/2, 5/2)}
           

SymPy 的集合論工具

集合論可以說是數學的基礎，在任何數學的方向上都能夠看到集合論的身影。在 SymPy 裡面，有一個類叫做 sympy.sets.sets.set。在集合論裡面，常見的就是邊界，補集，包含，并集，交集等常見的操作。但是感覺 SymPy 中的集合論操作主要集中在實數域或者複數域上。

對于閉區間

和開區間

而言，在 SymPy 中使用以下方法來表示：

I = sympy.Interval(0,1)
J = sympy.Interval.open(0,1)
K = sympy.Interval(0.5,2)
           

其開始和結束的點可以分别使用 start 和 end 來表示：

>>> I.start
0
>>> I.end
1
           

其長度用 measure 來表示：

>>> I.measure
1
           

其閉包用 closure 來表示：

>>> I.closure
Interval(0, 1)
           

其内點用 interior 來表示：

>>> I.interior
Interval.open(0, 1)
           

判斷其邊界條件可以使用 left_open 或者 right_open 來做：

>>> I.left_open
False
>>> I.right_open
False
           

對于兩個集合之間的操作，可以參考以下方法：

I = sympy.Interval(0,1)
K = sympy.Interval(0.5,2)
>>> I.intersect(K)
Interval(0.500000000000000, 1)
>>> I.union(K)
Interval(0, 2)
>>> I-K
Interval.Ropen(0, 0.500000000000000)
>>> K-I
Interval.Lopen(1, 2)
>>> I.symmetric_difference(K)
Union(Interval.Ropen(0, 0.500000000000000), Interval.Lopen(1, 2))
           

實數集

在 SymPy 中用 sympy.S.Reals 來表示，自然數使用 sympy.S.Naturals，非負整數用 sympy.S.Naturals0，整數用 sympy.S.Integers 來表示。補集的計算可以用減号，也可以使用 complement 函數。

>>> sympy.S.Reals
Reals
>>> sympy.S.Reals-I
Union(Interval.open(-oo, 0), Interval.open(1, oo))
>>> I.complement(sympy.S.Reals)
Union(Interval.open(-oo, 0), Interval.open(1, oo))
>>> sympy.S.Reals.complement(I)
EmptySet()
>>> I.complement(K)
Interval.Lopen(1, 2)
>>> I.complement(sympy.S.Reals)
Union(Interval.open(-oo, 0), Interval.open(1, oo))
           

SymPy 的邏輯工具

在邏輯運算中，我們可以使用

來代表元素。&, |, ~, >> 分别表示 AND，OR，NOT，imply。而邏輯運算同樣可以使用 sympy.simplify_logic 簡化。

A, B, C = sympy.symbols("A B C")
>>> sympy.simplify_logic(A | (A & B))
A
>>> sympy.simplify_logic((A>>B) & (B>>A))
(A & B) | (~A & ~B)
>>> A>>B
Implies(A, B)
           

SymPy 的級數工具

SymPy 的級數工具有一部分放在具體數學（Concrete Mathematics）章節了。有的時候，我們希望計算某個級數是發散的，還是收斂的，就可以使用 is_convergence 函數。考慮最常見的級數：

，

。

>>> n = sympy.Symbol("n", integer=True)
>>> sympy.Sum(1/n, (n,1,sympy.oo)).is_convergent()
False
>>> sympy.Sum(1/n**2, (n,1,sympy.oo)).is_convergent()
True
           

如果想計算出收斂級數的值，加上 doit() 函數即可；如果想計算有效數字，加上 evalf() 函數即可。

>>> sympy.Sum(1/n**2, (n,1,sympy.oo)).evalf()
1.64493406684823
>>> sympy.Sum(1/n**2, (n,1,sympy.oo)).doit()
pi**2/6
>>> sympy.Sum(1/n**3, (n,1,sympy.oo)).evalf()
1.20205690315959
>>> sympy.Sum(1/n**3, (n,1,sympy.oo)).doit()
zeta(3)
           

除了加法之外，SymPy 也支援連乘，其符号是 sympy.Product，考慮

，

。

>>> sympy.Product(n/(n+1), (n,1,sympy.oo)).is_convergent()
False
>>> sympy.Product(sympy.cos(sympy.pi/n), (n, 1, sympy.oo)).is_convergent()
True
           

SymPy 的 ODE 工具

在常微分方程（Ordinary Differential Equation）中，最常見的就是解方程，而解方程主要是靠 dsolve 函數。例如想求解以下的常微分方程：

，

，

，

可以使用 dsolve 函數：

>>> f = sympy.Function('f')
>>> sympy.dsolve(sympy.Derivative(f(x),x) + f(x), f(x))
Eq(f(x), C1*exp(-x))
>>> sympy.dsolve(sympy.Derivative(f(x),x,2) + f(x), f(x))
Eq(f(x), C1*sin(x) + C2*cos(x))
>>> sympy.dsolve(sympy.Derivative(f(x),x,3) + f(x), f(x))
Eq(f(x), C3*exp(-x) + (C1*sin(sqrt(3)*x/2) + C2*cos(sqrt(3)*x/2))*sqrt(exp(x)))
           

而常微分方程對于不同的方程類型也有着不同的解法，可以使用 classify_ode 來判斷常微分方程的類型：

>>> sympy.classify_ode(sympy.Derivative(f(x),x) + f(x), f(x))
('separable', '1st_exact', '1st_linear', 'almost_linear', '1st_power_series', 'lie_group', 'nth_linear_constant_coeff_homogeneous', 'separable_Integral', '1st_exact_Integral', '1st_linear_Integral', 'almost_linear_Integral')
>>> sympy.classify_ode(sympy.Derivative(f(x),x,2) + f(x), f(x))
('nth_linear_constant_coeff_homogeneous', '2nd_power_series_ordinary')
>>> sympy.classify_ode(sympy.Derivative(f(x),x,3) + f(x), f(x))
('nth_linear_constant_coeff_homogeneous',)
           

SymPy 的 PDE 工具

在偏微分方程（Partitial Differential Equation）中，同樣可以直接求解和判斷偏微分方程的類型，分别使用函數 pdsolve() 和 classify_pde()。假設

是一個二進制函數，分别滿足以下偏微分方程：

，

，

。

>>> f = sympy.Function("f")(x,y)
>>> sympy.pdsolve(sympy.Derivative(f,x)+sympy.Derivative(f,y),f)
Eq(f(x, y), F(x - y))
>>> sympy.pdsolve(f.diff(x)+f.diff(y)+f,f)
Eq(f(x, y), F(x - y)*exp(-x/2 - y/2))
>>> sympy.pdsolve(f.diff(x)+f.diff(y)+f+10,f)
Eq(f(x, y), F(x - y)*exp(-x/2 - y/2) - 10)
           

檢視類型就用 classify_pde() 函數：

>>> sympy.classify_pde(f.diff(x)+f.diff(y)+f)
('1st_linear_constant_coeff_homogeneous',)
>>> sympy.classify_pde(f.diff(x)+f.diff(y)+f+10,f)
('1st_linear_constant_coeff', '1st_linear_constant_coeff_Integral')
>>> sympy.classify_pde(f.diff(x)+f.diff(y)+f+10,f)
('1st_linear_constant_coeff', '1st_linear_constant_coeff_Integral')
           

不過目前的 PDE 解法貌似隻支援一階偏導數，二階或者以上的偏導數就不支援了。

SymPy 的數論工具

在數論中，素數就是一個最基本的概念之一。而素數的批量計算，比較快的方法就是篩法（sieve method）。在 sympy 中，同樣有 sympy.sieve 這個工具，用于計算素數。如果想輸出前100個素數，那麼

>>> sympy.sieve._reset()
>>> sympy.sieve.extend_to_no(100)
>>> sympy.sieve._list
array('l', [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97, 101, 103, 107, 109, 113, 127, 131, 137, 139, 149, 151, 157, 163, 167, 173, 179, 181, 191, 193, 197, 199, 211, 223, 227, 229, 233, 239, 241, 251, 257, 263, 269, 271, 277, 281, 283, 293, 307, 311, 313, 317, 331, 337, 347, 349, 353, 359, 367, 373, 379, 383, 389, 397, 401, 409, 419, 421, 431, 433, 439, 443, 449, 457, 461, 463, 467, 479, 487, 491, 499, 503, 509, 521, 523, 541, 547, 557, 563, 569, 571, 577, 587, 593, 599, 601, 607, 613, 617, 619, 631])
           

如果想輸出一個區間内的所有素數，可以使用 primerange(a,b) 函數：

>>> [i for i in sympy.sieve.primerange(10,100)]
[11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97]
           

search() 函數是為了計算某個數附近是第幾個素數：

>>> sympy.sieve.search(10)
(4, 5)
>>> sympy.sieve.search(11)
(5, 5)
           

如果隻想獲得第 n 個素數，則使用函數 sympy.ntheory.generate.prime(n) 即可。如果是希望計算 x 後面的下一個素數，使用 sympy.ntheory.generate.nextprime(x) 即可。判斷 x 是否是素數，可以使用 sympy.ntheory.generate.isprime(x)。

>>> sympy.ntheory.generate.prime(10)
29
>>> sympy.ntheory.generate.nextprime(10)
11
>>> sympy.ntheory.generate.nextprime(11)
13
>>> sympy.ntheory.generate.isprime(11)
True
>>> sympy.ntheory.generate.isprime(12)
False
           

除此之外，SymPy 的數論方法還有很多，需要讀者根據 SymPy 的官方文檔自行探索。

SymPy 的範疇論工具

SymPy 還支援範疇論（Category Theory）的一些計算方法，在這裡簡要地列舉一下。

>>> A = sympy.categories.Object("A")
>>> B = sympy.categories.Object("B")
>>> f = sympy.categories.NamedMorphism(A,B,"f")
>>> f.domain
Object("A")
>>> f.codomain
Object("B")
           

由于範疇論是數學的“黑話”，是以其餘方法留給範疇論的科研人員自行翻閱。

總結：

整體來看，

SymPy

是一個非常卓越的 Python 開源符号計算庫。在符号計算領域，不僅支援常見的微積分，線性代數，幾何運算，還支援集合論，微分方程，數論等諸多數學方向。後續筆者将會持續跟進并研究這一卓越的開源工具庫。

參考文獻：

1. Meurer A, Smith C P, Paprocki M, et al. SymPy: symbolic computing in Python[J]. PeerJ Computer Science, 2017, 3: e103.
2. GitHub：https://github.com/sympy/sympy
3. SymPy：https://www.sympy.org/en/index.html
4. Sympy 維基百科：https://en.wikipedia.org/wiki/SymPy
5. GreatX's Blog：數值 Python：符号計算：https://vlight.me/2018/04/01/Numerical-Python-Symbolic-Computing/
6. SymPy 符号計算-讓Python幫我們推公式：https://zhuanlan.zhihu.com/p/83822118
7. Python 科學計算利器---SymPy庫：https://www.jianshu.com/p/339c91ae9f41