Может кто-нибудь объяснить, почему scipy.integrate.quad дает разные результаты для одинаково больших диапазонов при интеграции sin (X)? - PullRequest
5 голосов
/ 24 февраля 2009

Я пытаюсь численно интегрировать произвольную (известную, когда я кодирую) функцию в моей программе используя методы численного интегрирования. Я использую 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)

Я считаю это поведение странным, потому что -
1. При обычной интеграции интеграция по полному циклу дает ноль.
2. При численном интегрировании это (1) не обязательно имеет место, потому что вы можете просто аппроксимируя общую площадь под кривой.

В любом случае, если предположить, что 1 - «Истина», или «2» - «Истина», я считаю, что такое поведение противоречиво. Либо обе интеграции (от -pi до pi и от 0 до 2 * pi) должны возвращать 0.0 (первое значение в кортеже - результат, а второе - ошибка), либо возвращать 2.257 ...

Может кто-нибудь объяснить, почему это происходит? Это действительно несоответствие? Может ли кто-то также сказать мне, если я что-то упускаю из-за численных методов?

В любом случае, в моем последнем приложении я планирую использовать вышеупомянутый метод, чтобы найти длину дуги функции. Если у кого-то есть опыт в этой области, пожалуйста, сообщите мне, как лучше всего это сделать на Python.

Редактировать
Примечание
У меня уже есть первые дифференциальные значения во всех точках диапазона, хранящихся в массиве.
Текущая ошибка допустима.
Конечная нота

Я читал об этом википедию. Как отметил Дмитрий, я буду интегрировать sqrt (1 + diff (f (x), x) ^ 2), чтобы получить длину дуги. Я хотел спросить: есть ли лучшее приближение / лучшие практики (?) / Более быстрый способ сделать это? Если потребуется больше контекста, я выложу его отдельно / опубликую контекст здесь, как вы пожелаете.

Ответы [ 6 ]

10 голосов
/ 24 февраля 2009

Функция quad - это функция из старой библиотеки Фортрана. Он работает, судя по плоскостности и наклону функции, которую он интегрирует, как обрабатывать размер шага, который он использует для численной интеграции, чтобы максимизировать эффективность. Это означает, что вы можете получить немного разные ответы от одного региона к другому, даже если они аналитически одинаковы.

Без сомнения, обе интеграции должны возвращать ноль. Возвращение чего-то 1 / (10 триллионов) довольно близко к нулю! Небольшие различия связаны с тем, что quad переворачивает sin и изменяет размеры шагов. Для вашего запланированного задания quad будет всем, что вам нужно.

EDIT: Для того, что вы делаете, я думаю, quad хорошо. Это быстро и довольно точно. Мое последнее утверждение - используйте его с уверенностью, если только вы не найдете то, что действительно пошло не так. Если он не возвращает бессмысленный ответ, то, вероятно, он работает просто отлично. Не беспокойся.

6 голосов
/ 24 февраля 2009

Я бы сказал, что число O (10 ^ -14) фактически равно нулю. Какова ваша терпимость?

Возможно, алгоритм, лежащий в основе четырехугольника, не самый лучший. Вы можете попробовать другой метод интеграции и посмотреть, улучшит ли это что-то. Рунге-Кутта 5-го порядка может быть очень хорошей техникой общего назначения.

Это может быть просто природа чисел с плавающей запятой: «Что должен знать каждый компьютерный специалист об арифметике с плавающей запятой».

6 голосов
/ 24 февраля 2009

Я думаю, что это, вероятно, точность машины, поскольку оба ответа фактически равны нулю.

Если вы хотите получить ответ изо рта лошади, я бы разместил этот вопрос на scipy форуме

4 голосов
/ 24 февраля 2009

Этот вывод мне кажется правильным, так как здесь у вас есть абсолютная оценка ошибки. Интегральное значение sin (x) действительно должно иметь значение ноль для полного периода (любой интервал 2 * pi длины) как в обычной, так и в числовой интеграции, и ваши результаты близки к этому значению.
Чтобы оценить длину дуги, вы должны вычислить интеграл для функции sqrt (1 + diff (f (x), x) ^ 2), где diff (f (x), x) является производной от f (x). Смотри также Длина дуги

3 голосов
/ 25 февраля 2009
0.0 == 2.3e-16 (absolute error tolerance 4.4e-14)

Оба ответа одинаковы и правильные , т. Е. Ноль в пределах данного допуска.

2 голосов
/ 04 декабря 2009

Разница проистекает из того факта, что sin (x) = - sin (-x) точно даже с конечной точностью. Принимая во внимание, что конечная точность только дает sin (x) ~ sin (x + 2 * pi) приблизительно. Конечно, было бы неплохо, если бы quad был достаточно умен, чтобы понять это, но на самом деле нет никакого способа узнать, априори, что интеграл по двум заданным вами интервалам эквивалентен или что первый результат лучше.

...