来自 Wester 文章的例子¶
介绍¶
在本教程中，我们展示了来自Wester文章中关于几个计算机代数系统数学能力比较和批评的例子（参见[Wester1999]_）。所有例子都与多项式和代数计算相关，并且为所有例子添加了SymPy特定的注释。
示例¶
本教程中的所有示例都是可计算的，因此可以直接将它们复制并粘贴到 Python shell 中，并利用它们做些有用的事情。所有计算都是使用以下设置完成的:
>>> from sympy import *
>>> init_printing(use_unicode=True)
>>> var('x,y,z,s,c')
(x, y, z, s, c)
简单的单变量多项式因式分解¶
要获得一个多项式的因式分解，请使用 factor() 函数。默认情况下，factor() 返回未计算的形式，因此输入多项式的内容保持未展开状态，如下例所示:
>>> factor(6*x - 10)
2⋅(3⋅x - 5)
要以更系统的方式达到相同的效果，请使用 primitive() 函数，该函数返回输入多项式的内容和原始部分:
>>> primitive(6*x - 10)
(2, 3⋅x - 5)
备注
内容和原始部分只能在一个环上计算。为了简化一个域上的多项式的系数，请使用 monic()。
单变量最大公约数、结式和因式分解¶
考虑整数上的单变量多项式 f、g 和 h:
>>> f = 64*x**34 - 21*x**47 - 126*x**8 - 46*x**5 - 16*x**60 - 81
>>> g = 72*x**60 - 25*x**25 - 19*x**23 - 22*x**39 - 83*x**52 + 54*x**10 + 81
>>> h = 34*x**19 - 25*x**16 + 70*x**7 + 20*x**3 - 91*x - 86
我们可以使用 gcd() 函数计算两个多项式的最大公约数（GCD）:
>>> gcd(f, g)
1
我们看到 f 和 g 没有公因数。然而，f*h 和 g*h 有一个明显的因数 h:
>>> gcd(expand(f*h), expand(g*h)) - h
0
同样可以使用一元多项式的结果来验证:
>>> resultant(expand(f*h), expand(g*h))
0
在整数上分解大的一元多项式（在这种情况下为120次）也是可能的:
>>> factor(expand(f*g))
⎛ 60 47 34 8 5 ⎞ ⎛ 60 52 39 25 23 10 ⎞
-⎝16⋅x + 21⋅x - 64⋅x + 126⋅x + 46⋅x + 81⎠⋅⎝72⋅x - 83⋅x - 22⋅x - 25⋅x - 19⋅x + 54⋅x + 81⎠
多元GCD和因式分解¶
在单变量情况下可以做的事情，在多变量多项式中也可以做。考虑以下在 \(\mathbb{Z}[x,y,z]\) 中的多项式 f、g 和 h:
>>> f = 24*x*y**19*z**8 - 47*x**17*y**5*z**8 + 6*x**15*y**9*z**2 - 3*x**22 + 5
>>> g = 34*x**5*y**8*z**13 + 20*x**7*y**7*z**7 + 12*x**9*y**16*z**4 + 80*y**14*z
>>> h = 11*x**12*y**7*z**13 - 23*x**2*y**8*z**10 + 47*x**17*y**5*z**8
如前所述，我们可以验证 f 和 g 没有公因数:
>>> gcd(f, g)
1
然而，f*h 和 g*h 有一个明显的因子 h:
>>> gcd(expand(f*h), expand(g*h)) - h
0
大型多项式的多元分解也是可能的:
>>> factor(expand(f*g))
7 ⎛ 9 9 3 7 6 5 12 7⎞ ⎛ 22 17 5 8 15 9 2 19 8 ⎞
-2⋅y ⋅z⋅⎝6⋅x ⋅y ⋅z + 10⋅x ⋅z + 17⋅x ⋅y⋅z + 40⋅y ⎠⋅⎝3⋅x + 47⋅x ⋅y ⋅z - 6⋅x ⋅y ⋅z - 24⋅x⋅y ⋅z - 5⎠
指数中的符号支持¶
由 sympy.polys 提供的多项式操作函数主要用于整数指数。然而，使用符号指数进行计算是完全有效的，例如:
>>> n = var('n')
>>> gcd(x**n - x**(2*n), x**n)
n
x
Results may depend on powers being expanded (which will depend on
assumptions of the base):
>>> gcd(x**(n + 4), x**(n + 1) + 3*x**n)
1
>>> x = var('x', positive=True)
>>> gcd(x**(n + 4), x**(n + 1) + 3*x**n)
n
x
测试多项式是否有公共零点¶
要测试两个多项式是否有共同的根，我们可以使用 resultant() 函数。理论表明，如果两个多项式有一个共同的零点，那么它们的结果为零。例如:
>>> x = var('x')
>>> resultant(3*x**4 + 3*x**3 + x**2 - x - 2, x**3 - 3*x**2 + x + 5)
0
我们可以通过分解多项式来可视化这一事实:
>>> factor(3*x**4 + 3*x**3 + x**2 - x - 2)
⎛ 3 ⎞
(x + 1)⋅⎝3⋅x + x - 2⎠
>>> factor(x**3 - 3*x**2 + x + 5)
⎛ 2 ⎞
(x + 1)⋅⎝x - 4⋅x + 5⎠
在这两种情况下，我们都得到了因子 \(x + 1\)，这告诉我们公共根是 \(x = -1\)。
归一化简单有理函数¶
要从有理函数的分子和分母中优雅地去除公因子，请使用 cancel() 函数。例如:
>>> cancel((x**2 - 4)/(x**2 + 4*x + 4))
x - 2
─────
x + 2
展开表达式和重新因式分解¶
可以在扩展形式和因式分解形式中轻松处理表达式。考虑一个以扩展形式表示的多项式 f。我们对它进行微分并将结果重新因式分解:
>>> f = expand((x + 1)**20)
>>> g = diff(f, x)
>>> factor(g)
19
20⋅(x + 1)
同样的效果可以通过因式分解形式实现:
>>> diff((x + 1)**20, x)
19
20⋅(x + 1)
以分圆多项式为项进行因式分解¶
SymPy 可以非常高效地分解形如 \(x^n \pm 1\) 的多项式，用循环多项式表示:
>>> factor(x**15 - 1)
⎛ 2 ⎞ ⎛ 4 3 2 ⎞ ⎛ 8 7 5 4 3 ⎞
(x - 1)⋅⎝x + x + 1⎠⋅⎝x + x + x + x + 1⎠⋅⎝x - x + x - x + x - x + 1⎠
原始的 Wester 例子是 \(x^{100} - 1\)，但为了可读性被截断了。注意，这对于 factor() 分解度为 1000 或更高的多项式来说并不困难。
高斯数上的单变量因式分解¶
考虑一个具有整数系数的一元多项式 f:
>>> f = 4*x**4 + 8*x**3 + 77*x**2 + 18*x + 153
我们希望在复数域上对 f 进行因式分解。为此，我们像之前一样使用 factor()，但这次我们将 gaussian 关键字设置为 True:
>>> factor(f, gaussian=True)
⎛ 3⋅ⅈ⎞ ⎛ 3⋅ⅈ⎞
4⋅⎜x - ───⎟⋅⎜x + ───⎟⋅(x + 1 - 4⋅ⅈ)⋅(x + 1 + 4⋅ⅈ)
⎝ 2 ⎠ ⎝ 2 ⎠
结果我们得到了一个 f 的因式分解，其因子是首一多项式（这是在 SymPy 中计算域时的一般规则）。gaussian 关键字有助于提高代码的可读性，但同样的结果可以使用更一般的语法计算:
>>> factor(f, extension=I)
⎛ 3⋅ⅈ⎞ ⎛ 3⋅ⅈ⎞
4⋅⎜x - ───⎟⋅⎜x + ───⎟⋅(x + 1 - 4⋅ⅈ)⋅(x + 1 + 4⋅ⅈ)
⎝ 2 ⎠ ⎝ 2 ⎠
使用自动字段扩展进行计算¶
考虑两个单变量多项式 f 和 g:
>>> f = x**3 + (sqrt(2) - 2)*x**2 - (2*sqrt(2) + 3)*x - 3*sqrt(2)
>>> g = x**2 - 2
我们希望减少有理函数 f/g 的分子和分母的次数。为此，我们使用 cancel() 函数:
>>> cancel(f/g)
3 2 2
x - 2⋅x + √2⋅x - 3⋅x - 2⋅√2⋅x - 3⋅√2
───────────────────────────────────────
2
x - 2
不幸的是，没有发生任何有趣的事情。这是因为默认情况下，SymPy 将 \(\sqrt{2}\) 视为一个生成器，从而得到分子上的一个二元多项式。要让 cancel() 识别 \(\sqrt{2}\) 的代数性质，需要使用 extension 关键字:
>>> cancel(f/g, extension=True)
2
x - 2⋅x - 3
────────────
x - √2
设置 extension=True 告诉 cancel() 为 f/g 的系数找到最小的代数数域。自动推断的域是 \(\mathbb{Q}(\sqrt{2})\)。如果不希望依赖自动推断，可以通过设置 extension 关键字并指定一个明确的代数数来获得相同的结果:
>>> cancel(f/g, extension=sqrt(2))
2
x - 2⋅x - 3
────────────
x - √2
在各种域上的单变量因式分解¶
考虑一个具有整数系数的一元多项式 f:
>>> f = x**4 - 3*x**2 + 1
使用 sympy.polys 我们可以获得 f 在不同域上的因式分解，其中包括：
rationals:
>>> factor(f) ⎛ 2 ⎞ ⎛ 2 ⎞ ⎝x - x - 1⎠⋅⎝x + x - 1⎠
有限域:
>>> factor(f, modulus=5) 2 2 (x - 2) ⋅(x + 2)
代数数:
>>> alg = AlgebraicNumber((sqrt(5) - 1)/2, alias='alpha') >>> factor(f, extension=alg) (x - α)⋅(x + α)⋅(x - 1 - α)⋅(x + α + 1)
将多项式分解为线性因子¶
目前，SymPy 可以将多项式分解为各种域上的不可约多项式，这可能导致分裂分解（分解为线性因子）。然而，目前没有系统的方法来自动推断分裂域（代数数域）。未来将实现以下语法:
>>> factor(x**3 + x**2 - 7, split=True)
Traceback (most recent call last):
...
NotImplementedError: 'split' option is not implemented yet
注意这与 extension=True 不同，因为后者仅指示表达式解析应如何进行，而不是计算域应是什么。可以使用 solve() 函数为几类多项式模拟 split 关键字。
有限域上的高级因式分解¶
考虑一个具有整数系数的一元多项式 f:
>>> f = x**11 + x + 1
我们可以在一个大的有限域 \(F_{65537}\) 上对 f 进行因式分解:
>>> factor(f, modulus=65537)
⎛ 2 ⎞ ⎛ 9 8 6 5 3 2 ⎞
⎝x + x + 1⎠⋅⎝x - x + x - x + x - x + 1⎠
并将生成的因式分解展开回:
>>> expand(_)
11
x + x + 1
获取多项式 f。这是通过有限域上的对称多项式表示完成的。同样的事情可以通过非对称表示来完成:
>>> factor(f, modulus=65537, symmetric=False)
⎛ 2 ⎞ ⎛ 9 8 6 5 3 2 ⎞
⎝x + x + 1⎠⋅⎝x + 65536⋅x + x + 65536⋅x + x + 65536⋅x + 1⎠
与对称表示一样，我们可以扩展因式分解以得到输入的多项式。然而，这次我们需要将扩展多项式的系数截断模65537:
>>> trunc(expand(_), 65537)
11
x + x + 1
使用多项式表达式¶
考虑一个多元多项式 f 在 \(\mathbb{Z}[x,y,z]\):
>>> f = expand((x - 2*y**2 + 3*z**3)**20)
我们想要计算 f 的因式分解。为此，我们通常使用 factor，但需要注意的是，考虑的多项式已经是展开形式，因此我们可以告诉因式分解例程跳过展开 f:
>>> factor(f, expand=False)
20
⎛ 2 3⎞
⎝x - 2⋅y + 3⋅z ⎠
在 sympy.polys 中，默认情况下会将所有作为参数传递给多项式操作函数和 Poly 类的表达式展开。如果我们知道展开是不必要的，那么通过设置 expand=False，我们可以为复杂的输入节省相当多的时间。这在处理以下表达式时尤为重要:
>>> g = expand((sin(x) - 2*cos(y)**2 + 3*tan(z)**3)**20)
>>> factor(g, expand=False)
20
⎛ 2 3 ⎞
⎝-sin(x) + 2⋅cos (y) - 3⋅tan (z)⎠
计算约化Gröbner基¶
要计算一组多项式的简化Gröbner基，请使用 groebner() 函数。该函数接受各种单项式排序，例如：lex、grlex 和 grevlex，或通过 order 关键字定义的用户自定义排序。lex 排序是最有趣的，因为它具有消除性质，这意味着如果传递给 groebner() 的多项式方程系统是零维的（有有限个解），则基的最后一个元素是一个单变量多项式。考虑以下示例:
>>> f = expand((1 - c**2)**5 * (1 - s**2)**5 * (c**2 + s**2)**10)
>>> groebner([f, c**2 + s**2 - 1])
⎛⎡ 2 2 20 18 16 14 12 10⎤ ⎞
GroebnerBasis⎝⎣c + s - 1, c - 5⋅c + 10⋅c - 10⋅c + 5⋅c - c ⎦, s, c, domain=ℤ, order=lex⎠
结果是一个普通的 Python 列表，因此我们可以轻松地对其所有元素应用一个函数，例如我们可以对这些元素进行因式分解:
>>> list(map(factor, _))
⎡ 2 2 10 5 5⎤
⎣c + s - 1, c ⋅(c - 1) ⋅(c + 1) ⎦
从上述内容中，我们可以轻松找到多项式方程组的所有解。或者，我们可以使用 solve() 以更系统的方式实现这一点:
>>> solve([f, s**2 + c**2 - 1], c, s)
[(-1, 0), (0, -1), (0, 1), (1, 0)]
代数数上的多元因式分解¶
在各种域上计算多元多项式与单变量情况一样简单。例如，考虑以下在 \(\mathbb{Q}(\sqrt{-3})\) 上的因式分解:
>>> factor(x**3 + y**3, extension=sqrt(-3))
⎛ ⎛ 1 √3⋅ⅈ⎞⎞ ⎛ ⎛ 1 √3⋅ⅈ⎞⎞
(x + y)⋅⎜x + y⋅⎜- ─ - ────⎟⎟⋅⎜x + y⋅⎜- ─ + ────⎟⎟
⎝ ⎝ 2 2 ⎠⎠ ⎝ ⎝ 2 2 ⎠⎠
备注
目前不支持有限域上的多元多项式。
部分分式分解¶
考虑一个具有整数系数的一元有理函数 f:
>>> f = (x**2 + 2*x + 3)/(x**3 + 4*x**2 + 5*x + 2)
要将 f 分解为部分分式，请使用 apart() 函数:
>>> apart(f)
3 2 2
───── - ───── + ────────
x + 2 x + 1 2
(x + 1)
要从部分分式返回到有理函数，请使用 together() 和 cancel() 的组合:
>>> cancel(together(_))
2
x + 2⋅x + 3
───────────────────
3 2
x + 4⋅x + 5⋅x + 2
文献¶
Michael J. Wester, 《对CA系统数学能力的批判》, 1999, https://www.math.unm.edu/~wester/cas/book/Wester.pdf