python - 使用 scipy.integrate.simps 进行集成

Question

我正在尝试了解 scipy 包，但遇到了一些我无法理解的东西。

from scipy.integrate import simps
import numpy as np
def f1(x):
    ...:     return x**2
x = np.array([1,3,4])
y1 = f1(x)
I1 = integrate.simps(y1,x)
print(I1)
21.0

这完全对应于

∫₁⁴ x² dx = 21，

我没有得到的是 x = np.array([1, 3, 4]) 行。为什么这里需要 3 个？ 1和4是积分的极限，那么3是多少呢？有人可以向我解释一下吗？

Answer 1

documentation scipy.integrate.simps 的 说:

y : array_like

    Array to be integrated.

x : array_like, optional

    If given, the points at which y is sampled.

这些是对函数进行采样的点。由于您没有将要集成的函数直接传递给算法，因此您必须提供样本点。第二个数组给出您在上一行中计算的 y 值 y1 的 x 位置。尽管数值积分方法的某些实现直接采用被积函数，但它们始终会创建像您在此处提供的样本点。

因此数组x不是积分区间，尽管它的最大值和最小值给出了区间。

一般来说，对于任何数值积分算法，分布在积分区间上的样本点数量较多，会提高数值结果的准确性，而只有 3 个点几乎肯定会导致非常差的近似值。

但是，在您的示例中，被积函数是一个简单的 2 阶多项式。这样很容易积分(无论是分析还是数值)。您在 scipy.integrate.simps 中使用的算法是 Simpson's rule ，它基于将被积函数扩展至 2 阶。因此该方法能够准确求解您的样本积分。

要完全定义二阶多项式，您需要指定 3 个系数，并且为了能够导出这些系数，算法需要知道二阶多项式的至少 3 个点。然而，附加的第四点不会提供任何更多信息，因为曲线已经由三个点完全指定。这就是为什么在此示例中 3 个点足以给出准确结果的原因。

如果您不提供带有示例位置的列表 x，结果通常会是错误的，因为 y1< 中各个 y 值之间的间距将被假定为 1/。 (请参阅上面的文档链接)

另外，作为旁注，使用 Python 2.7 与 Numpy 1.6 和 Scipy 0.10 ，上面代码的结果是 20.75，可能是因为 x 的类型被假定为整数。明确声明它们是 float 的

x = np.array([1.0,3.0,4.0])

解决了这个问题，结果始终是21.0。您还可以看到，只要在 1.0 和 4.0 之间，实际的中间值并不重要。