NumPy参考 >例行程序 >Polynomials >Polynomial Package > Using the Convenience Classes
多项式程序包提供的便利类为:
名称
提供
多项式
电源系列
切比雪夫
切比雪夫系列
勒让德
勒让德系列
拉盖尔
拉盖尔系列
赫米特人
埃尔米特系列
埃尔米特
HermiteE系列
在这种情况下,级数是相应的多项式基函数的有限总和乘以系数。例如,幂级数看起来像
并且具有系数。具有相同系数的契比雪夫系列看起来像
更普遍地
在这种情况下,这里是度的Chebyshev函数,但也很容易成为其他任何类的基础函数。所有类的约定是系数与度i的基函数匹配。
所有这些类都是不可变的,并且具有相同的方法,尤其是它们实现了Python数字运算符+,-,*,//,%,divmod,**,==和!=。由于浮点舍入错误,后两个可能会有些问题。现在,我们快速演示使用NumPy 1.7.0版进行的各种操作。
首先,我们需要一个多项式类和一个多项式实例来使用。这些类可以直接从多项式包或相关类型的模块中导入。在这里,由于其熟悉性,我们从包中导入并使用常规的Polynomial类:
>>> from numpy.polynomial import Polynomial as P
>>> p = P([1,2,3])
>>> p
Polynomial([ 1., 2., 3.], domain=[-1, 1], window=[-1, 1])
请注意,打印输出的长版包含三个部分。第一个是系数,第二个是域,第三个是窗口:
>>> p.coef
array([ 1., 2., 3.])
>>> p.domain
array([-1., 1.])
>>> p.window
array([-1., 1.])
打印多项式会产生一个没有域和窗口的较短形式:
>>> print p
poly([ 1. 2. 3.])
当我们进行拟合时,我们将处理域和窗口,此刻我们将忽略它们并进行基本的代数和算术运算。
加减:
>>> p + p
Polynomial([ 2., 4., 6.], domain=[-1, 1], window=[-1, 1])
>>> p - p
Polynomial([ 0.], domain=[-1, 1], window=[-1, 1])
乘法:
>>> p * p
Polynomial([ 1., 4., 10., 12., 9.], domain=[-1, 1], window=[-1, 1])
权力:
>>> p**2
Polynomial([ 1., 4., 10., 12., 9.], domain=[-1, 1], window=[-1, 1])
师:
地板除法“ //”是多项式类别的除法运算符,在这方面,多项式被视为整数。对于Python版本3.x以下版本,“ /”运算符映射为“ //”,就像Python一样,对于更高版本,“ /”仅适用于按标量除法。在某些时候,它会被弃用:
>>> p // P([-1, 1])
Polynomial([ 5., 3.], domain=[-1, 1], window=[-1, 1])
余:
>>> p % P([-1, 1])
Polynomial([ 6.], domain=[-1, 1], window=[-1, 1])
Divmod:
>>> quo, rem = divmod(p, P([-1, 1]))
>>> quo
Polynomial([ 5., 3.], domain=[-1, 1], window=[-1, 1])
>>> rem
Polynomial([ 6.], domain=[-1, 1], window=[-1, 1])
评价:
>>> x = np.arange(5)
>>> p(x)
array([ 1., 6., 17., 34., 57.])
>>> x = np.arange(6).reshape(3,2)
>>> p(x)
array([[ 1., 6.],
[ 17., 34.],
[ 57., 86.]])
代换:
将多项式代入x并扩展结果。在这里,我们将p本身替换为扩展后导致新的4级多项式的p。如果将多项式视为函数,则这是函数的组合:
>>> p(p)
Polynomial([ 6., 16., 36., 36., 27.], domain=[-1, 1], window=[-1, 1])
根:
>>> p.roots()
array([-0.33333333-0.47140452j, -0.33333333+0.47140452j])
显式使用多项式实例并不总是很方便,因此元组,列表,数组和标量会在算术运算中自动进行转换:
>>> p + [1, 2, 3]
Polynomial([ 2., 4., 6.], domain=[-1, 1], window=[-1, 1])
>>> [1, 2, 3] * p
Polynomial([ 1., 4., 10., 12., 9.], domain=[-1, 1], window=[-1, 1])
>>> p / 2
Polynomial([ 0.5, 1. , 1.5], domain=[-1, 1], window=[-1, 1])
域,窗口或类不同的多项式不能在算术中混合:
>>> from numpy.polynomial import Chebyshev as T
>>> p + P([1], domain=[0,1])
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "<string>", line 213, in __add__
TypeError: Domains differ
>>> p + P([1], window=[0,1])
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "<string>", line 215, in __add__
TypeError: Windows differ
>>> p + T([1])
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "<string>", line 211, in __add__
TypeError: Polynomial types differ
但是可以使用不同的类型进行替换。实际上,这是如何对类型,域和窗口转换进行多项式类之间的转换:
>>> p(T([0, 1]))
Chebyshev([ 2.5, 2. , 1.5], domain=[-1, 1], window=[-1, 1])
这给出了Chebyshev形式的多项式p。这样做是因为 ,而代以对不改变原多项式。但是,所有乘法和除法操作都将使用Chebyshev级数进行,因此将得出结果的类型。
意图是所有多项式实例都是不可变的,因此有意地未实现扩展操作(+=
,-=
等)以及违反多项式实例的不可变性的任何其他功能。
多项式实例可以集成和区分。
>>> from numpy.polynomial import Polynomial as P
>>> p = P([2, 6])
>>> p.integ()
Polynomial([ 0., 2., 3.], domain=[-1, 1], window=[-1, 1])
>>> p.integ(2)
Polynomial([ 0., 0., 1., 1.], domain=[-1, 1], window=[-1, 1])
第一个示例对p进行一次积分,第二个示例对p进行两次积分。缺省情况下,积分的下界和积分常数为0,但都可以指定。
>>> p.integ(lbnd=-1)
Polynomial([-1., 2., 3.], domain=[-1, 1], window=[-1, 1])
>>> p.integ(lbnd=-1, k=1)
Polynomial([ 0., 2., 3.], domain=[-1, 1], window=[-1, 1])
在第一种情况下,积分的下限设置为-1,积分常数为0。在第二种情况下,积分的常数也设置为1。微分比较简单,因为唯一的选择是多项式被微分的次数:
>>> p = P([1, 2, 3])
>>> p.deriv(1)
Polynomial([ 2., 6.], domain=[-1, 1], window=[-1, 1])
>>> p.deriv(2)
Polynomial([ 6.], domain=[-1, 1], window=[-1, 1])
通过指定系数构造多项式只是获得多项式实例的一种方法,也可以通过指定其根,通过从其他多项式类型进行转换以及通过最小二乘拟合来创建多项式。拟合将在其单独的部分中讨论,其他方法如下所示:
>>> from numpy.polynomial import Polynomial as P
>>> from numpy.polynomial import Chebyshev as T
>>> p = P.fromroots([1, 2, 3])
>>> p
Polynomial([ -6., 11., -6., 1.], domain=[-1, 1], window=[-1, 1])
>>> p.convert(kind=T)
Chebyshev([ -9. , 11.75, -3. , 0.25], domain=[-1, 1], window=[-1, 1])
convert方法还可以转换域和窗口:
>>> p.convert(kind=T, domain=[0, 1])
Chebyshev([-2.4375 , 2.96875, -0.5625 , 0.03125], [ 0., 1.], [-1., 1.])
>>> p.convert(kind=P, domain=[0, 1])
Polynomial([-1.875, 2.875, -1.125, 0.125], [ 0., 1.], [-1., 1.])
在numpy的版本> = 1.7.0的基础和铸造类的方法也是可用的。强制转换方法的工作方式类似于转换方法,而基本方法则返回给定度数的基本多项式:
>>> P.basis(3)
Polynomial([ 0., 0., 0., 1.], domain=[-1, 1], window=[-1, 1])
>>> T.cast(p)
Chebyshev([ -9. , 11.75, -3. , 0.25], domain=[-1, 1], window=[-1, 1])
类型之间的转换可能很有用,但不建议日常使用。从Chebyshev阶数为50的级数到相同阶数的多项式级数时,数值精度的损失会使数值评估的结果基本上是随机的。
拟合是domain和window属性属于便捷类的一部分的原因。为了说明这个问题,下面绘制了直至5级的Chebyshev多项式的值。
>>> import matplotlib.pyplot as plt
>>> from numpy.polynomial import Chebyshev as T
>>> x = np.linspace(-1, 1, 100)
>>> for i in range(6): ax = plt.plot(x, T.basis(i)(x), lw=2, label="$T_%d$"%i)
...
>>> plt.legend(loc="upper left")
<matplotlib.legend.Legend object at 0x3b3ee10>
>>> plt.show()
在-1 <= x <= 1 范围内,它们是不错的,等波动函数介于+/- 1之间。在-2 <= x <= 2 范围内的相同图看起来非常不同:
>>> import matplotlib.pyplot as plt
>>> from numpy.polynomial import Chebyshev as T
>>> x = np.linspace(-2, 2, 100)
>>> for i in range(6): ax = plt.plot(x, T.basis(i)(x), lw=2, label="$T_%d$"%i)
...
>>> plt.legend(loc="lower right")
<matplotlib.legend.Legend object at 0x3b3ee10>
>>> plt.show()
可以看出,“好”部分已经缩小到无关紧要的程度。在使用Chebyshev多项式进行拟合时,我们要使用x在-1和1之间的区域,这就是窗口指定的区域。但是,要拟合的数据不可能在该时间间隔内拥有所有数据点,因此我们使用域指定数据点所在的间隔。拟合完成后,首先通过线性变换将域映射到窗口,然后使用映射的数据点完成通常的最小二乘拟合。拟合的窗口和范围是返回系列的一部分,在计算值,导数等时会自动使用。如果未在调用中指定它们,则拟合例程将使用默认窗口和包含所有数据点的最小域。如下所示,以拟合噪声正弦曲线。
>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> from numpy.polynomial import Chebyshev as T
>>> np.random.seed(11)
>>> x = np.linspace(0, 2*np.pi, 20)
>>> y = np.sin(x) + np.random.normal(scale=.1, size=x.shape)
>>> p = T.fit(x, y, 5)
>>> plt.plot(x, y, 'o')
[<matplotlib.lines.Line2D object at 0x2136c10>]
>>> xx, yy = p.linspace()
>>> plt.plot(xx, yy, lw=2)
[<matplotlib.lines.Line2D object at 0x1cf2890>]
>>> p.domain
array([ 0. , 6.28318531])
>>> p.window
array([-1., 1.])
>>> plt.show()