python 有人可以解释为什么 scipy.integrate.quad 在积分 sin(X) 时对同样长的范围给出不同的结果吗？

Can someone explain why scipy.integrate.quad gives different results for equally long ranges while integrating sin(X)?

提问by batbrat

I am trying to numerically integrate an arbitrary (known when I code) function in my program using numerical integration methods. I am using Python 2.5.2 along with SciPy's numerical integration package. In order to get a feel for it, i decided to try integrating sin(x) and observed this behavior-

我正在尝试使用数值积分方法在我的程序中对任意（在我编码时已知）函数进行数值积分。我正在使用 Python 2.5.2 和 SciPy 的数值集成包。为了感受一下，我决定尝试积分 sin(x) 并观察到这种行为-

>>> from math import pi
>>> from scipy.integrate import quad
>>> from math import sin
>>> def integrand(x):
...     return sin(x)
... 
>>> quad(integrand, -pi, pi)
(0.0, 4.3998892617846002e-14)
>>> quad(integrand, 0, 2*pi)
(2.2579473462709165e-16, 4.3998892617846002e-14)

I find this behavior odd because -
1. In ordinary integration, integrating over the full cycle gives zero.
2. In numerical integration, this (1) isn't necessarily the case, because you may just be approximating the total area under the curve.

我觉得这种行为很奇怪，因为 -
1. 在普通积分中，在整个周期内积分为零。
2. 在数值积分中，这 (1) 不一定是这种情况，因为您可能只是在近似曲线下的总面积。

In any case, either assuming 1 is True or assuming 2 is True, I find the behavior to be inconsistent. Either both integrations (-pi to pi and 0 to 2*pi) should return 0.0 (first value in the tuple is the result and the second is the error) or return 2.257...

在任何情况下，假设 1 为 True 或假设 2 为 True，我发现行为不一致。两个积分（-pi 到 pi 和 0 到 2*pi）都应该返回 0.0（元组中的第一个值是结果，第二个值是错误）或返回 2.257...

Can someone please explain why this is happening? Is this really an inconsistency? Can someone also tell me if I am missing something really basic about numerical methods?

有人可以解释为什么会这样吗？这真的是矛盾吗？有人也可以告诉我我是否遗漏了一些关于数值方法的基本知识？

In any case, in my final application, I plan to use the above method to find the arc length of a function. If someone has experience in this area, please advise me on the best policy for doing this in Python.

无论如何，在我最终的应用程序中，我打算使用上述方法来求一个函数的弧长。如果有人在这方面有经验，请告诉我在 Python 中执行此操作的最佳策略。

Edit
Note
I already have the first differential values at all points in the range stored in an array.
Current error is tolerable.
End note

编辑
注意
我已经在数组中存储了范围内所有点的第一个差分值。
当前的错误是可以容忍的。
尾注

I have read Wikipaedia on this. As Dimitry has pointed out, I will be integrating sqrt(1+diff(f(x), x)^2) to get the Arc Length. What I wanted to ask was - is there a better approximation/ Best practice(?) / faster way to do this. If more context is needed, I'll post it separately/ post context here, as you wish.

我已经阅读了维基百科。正如 Dimitry 所指出的，我将整合 sqrt(1+diff(f(x), x)^2) 以获得弧长。我想问的是 - 有没有更好的近似/最佳实践（？）/更快的方法来做到这一点。如果需要更多上下文，我会根据您的意愿单独发布/在此处发布上下文。

回答by physicsmichael

The quadfunction is a function from an old Fortran library. It works by judging by the flatness and slope of the function it is integrating how to treat the step size it uses for numerical integration in order to maximize efficiency. What this means is that you may get slightly different answers from one region to the next even if they're analytically the same.

该quad函数是来自旧 Fortran 库的函数。它的工作原理是根据它正在积分的函数的平坦度和斜率来判断如何处理它用于数值积分的步长，以最大限度地提高效率。这意味着即使在分析上相同，您可能会从一个区域到另一个区域得到略有不同的答案。

Without a doubt both integrations should return zero. Returning something that is 1/(10 trillion) is pretty close to zero! The slight differences are due to the way quadis rolling over sinand changing its step sizes. For your planned task, quadwill be all you need.

毫无疑问，两个积分都应该返回零。返回 1/（10 万亿）的东西非常接近于零！细微的差异是由于方式quad正在滚动sin并改变其步长。对于您计划的任务，quad将是您所需要的。

EDIT: For what you're doing I think quadis fine. It is fast and pretty accurate. My final statement is use it with confidence unless you find something that really has gone quite awry. If it doesn't return a nonsensical answer then it is probably working just fine. No worries.

编辑：对于你在做什么，我认为quad很好。它很快而且非常准确。我的最后声明是自信地使用它，除非您发现某些东西确实出了问题。如果它没有返回一个荒谬的答案，那么它可能工作得很好。不用担心。

回答by Simon

I think it is probably machine precision since both answers are effectively zero.

我认为这可能是机器精度，因为两个答案实际上都为零。

If you want an answer from the horse's mouth I would post this question on the scipy discussion board

如果你想从马嘴里得到答案，我会在scipy 讨论板上发布这个问题

回答by duffymo

I would say that a number O(10^-14) is effectively zero. What's your tolerance?

我会说数字 O(10^-14) 实际上为零。你的容忍度是多少？

It might be that the algorithm underlying quad isn't the best. You might try another method for integration and see if that improves things. A 5th order Runge-Kutta can be a very nice general purpose technique.

可能是 quad 底层的算法不是最好的。您可以尝试另一种集成方法，看看是否能改善情况。五阶 Runge-Kutta 可以是一种非常好的通用技术。

It could be just the nature of floating point numbers: "What Every Computer Scientist Should Know About Floating Point Arithmetic".

这可能只是浮点数的本质：“每个计算机科学家都应该了解浮点运算”。

回答by okutane

This output seems correct to me since you have absolute error estimate here. The integral value of sin(x) is indeed should have value of zero for full period (any interval of 2*pi length) in both ordinary and numeric integration and your results is close to that value.
To evaluate arc length you should calculate integral for sqrt(1+diff(f(x), x)^2) function, where diff(f(x), x) is derivative of f(x). See also Arc length

这个输出对我来说似乎是正确的，因为你在这里有绝对误差估计。sin(x) 的积分值在普通积分和数字积分中的整个周期（任何 2*pi 长度的间隔）的值确实应该为零，并且您的结果接近该值。
要评估弧长，您应该计算 sqrt(1+diff(f(x), x)^2) 函数的积分，其中 diff(f(x), x) 是 f(x) 的导数。另见弧长

回答by jfs

0.0 == 2.3e-16 (absolute error tolerance 4.4e-14)

Both answers are the same and correcti.e., zero within the given tolerance.

两个答案相同且正确，即在给定容差内为零。

回答by Kevin Mitchell

The difference comes from the fact that sin(x)=-sin(-x) exactly even in finite precision. Whereas finite precision only gives sin(x)~sin(x+2*pi) approximately. Sure it would be nice if quad were smart enough to figure this out, but it really has no way of knowing apriori that the integral over the two intervals you give are equivalent or that the the first is a better result.

不同之处在于 sin(x)=-sin(-x) 即使在有限精度下也是如此。而有限精度仅给出 sin(x)~sin(x+2*pi) 近似值。当然，如果四边形足够聪明来解决这个问题会很好，但它确实无法先验地知道你给出的两个区间的积分是等价的，或者第一个是更好的结果。