math::calculus::romberg

ИМЯ

math::calculus::romberg — интегрирование методом Ромберга

Содержание

Содержание
Краткое описание
Описание
ПРОЦЕДУРЫ
ПАРАМЕТРЫ
ПАРАМЕТРЫ
ОПИСАНИЕ
НЕСОБСТВЕННЫЕ ИНТЕГРАЛЫ
ДРУГИЕ ЗАМЕНЫ ПЕРЕМЕННОЙ
Ошибки, идеи, отзывы
См. также
Категория
Авторские права

КРАТКОЕ ОПИСАНИЕ

package require Tcl 8.5 9
package require math::calculus 0.6

ОПИСАНИЕ

Процедуры romberg пакета math::calculus выполняют численное интегрирование функции одной переменной. Они предназначены для использования в реальных приложениях: отличаются надёжностью и точностью, а также обеспечивают приемлемую эффективность по количеству вычислений функции.

ПРОЦЕДУРЫ

Для интегрирования методом Ромберга доступны следующие процедуры:

::math::calculus::romberg f a b ?-option value...?

Вычисляет интеграл аналитической функции на заданном интервале.
::math::calculus::romberg_infinity f a b ?-option value...?

Вычисляет интеграл аналитической функции на полунеограниченном интервале.
::math::calculus::romberg_sqrtSingLower f a b ?-option value...?

Вычисляет интеграл функции, которая, как предполагается, аналитична на интервале, за исключением обратной корневой особенности в нижнем пределе.
::math::calculus::romberg_sqrtSingUpper f a b ?-option value...?

Вычисляет интеграл функции, которая, как предполагается, аналитична на интервале, за исключением обратной корневой особенности в верхнем пределе.
::math::calculus::romberg_powerLawLower gamma f a b ?-option value...?

Вычисляет интеграл функции, которая, как предполагается, аналитична на интервале, за исключением степенной особенности в нижнем пределе.
::math::calculus::romberg_powerLawUpper gamma f a b ?-option value...?

Вычисляет интеграл функции, которая, как предполагается, аналитична на интервале, за исключением степенной особенности в верхнем пределе.
::math::calculus::romberg_expLower f a b ?-option value...?

Вычисляет интеграл экспоненциально возрастающей функции; нижний предел области интегрирования может быть сколь угодно большим отрицательным числом.
::math::calculus::romberg_expUpper f a b ?-option value...?

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

ПАРАМЕТРЫ

f

Функция для интегрирования. Должна быть задана одной командой Tcl, к которой будет добавлен один аргумент — абсцисса, в которой вычисляется функция. Перед вычислением первое слово команды обрабатывается с помощью namespace which в контексте вызывающей процедуры. Благодаря этой обработке команда может быть локальной для вызывающего пространства имён, а не обязательно глобальной.
a

Нижний предел области интегрирования.
b

Верхний предел области интегрирования. Для процедур romberg_sqrtSingLower, romberg_sqrtSingUpper, romberg_powerLawLower, romberg_powerLawUpper, romberg_expLower и romberg_expUpper нижний предел должен быть строго меньше верхнего. Для остальных процедур пределы могут быть заданы в любом порядке.
gamma

Показатель степени для степенной особенности; подробности см. в разделе НЕСОБСТВЕННЫЕ ИНТЕГРАЛЫ.

ПАРАМЕТРЫ

-abserror epsilon

Задаёт условие остановки интегрирования: оно продолжается, пока оценка абсолютной погрешности интеграла не станет меньше epsilon. Если функция (или любая из её производных) имеет особенности, погрешность может быть значительно завышена или занижена; подробности см. в разделе НЕСОБСТВЕННЫЕ ИНТЕГРАЛЫ. Значение по умолчанию — 1.0e-08.
-relerror epsilon

Задаёт условие остановки интегрирования: оно продолжается, пока оценка относительной погрешности интеграла не станет меньше epsilon. Если функция (или любая из её производных) имеет особенности, погрешность может быть значительно завышена или занижена; подробности см. в разделе НЕСОБСТВЕННЫЕ ИНТЕГРАЛЫ. Значение по умолчанию — 1.0e-06.
-maxiter m

Задаёт условие завершения интегрирования после не более чем n троекратных увеличений количества выполненных вычислений. Иными словами, при значении n для -maxiter механизм интегрирования выполнит не более 3**n вычислений функции. Значение по умолчанию — 14, что соответствует пределу примерно в 4,8 миллиона вычислений. (Для хорошо ведущих себя функций обычно достаточно нескольких сотен вычислений.)
-degree d

Задаёт степень d экстраполирующего полинома, используемого при интегрировании методом Ромберга; подробности см. в разделе ОПИСАНИЕ. Значение по умолчанию — 4. Максимальное значение — m-1.

ОПИСАНИЕ

Процедура romberg выполняет интегрирование методом Ромберга с использованием модифицированного правила средней точки. Интегрирование методом Ромберга — итерационный процесс. На первом шаге функция вычисляется в середине области интегрирования, а полученное значение умножается на ширину интервала, что даёт наиболее грубую оценку. На втором шаге интервал делится на три части, функция вычисляется в середине каждой части, а сумма значений умножается на три. На третьем шаге используются девять частей, на четвёртом — двадцать семь и так далее: на каждом шаге количество подынтервалов увеличивается втрое.

После того как интервал разделён не менее чем d раз, по интегралам, оценённым при последних d+1 разбиениях, строится полином. Считается, что интегралы являются функцией квадрата ширины подынтервалов (этот процесс описывается в любом хорошем учебнике по численному анализу в разделе «Интегрирование методом Ромберга»). Полином экстраполируется к нулевому шагу, чтобы вычислить значение интеграла и оценку погрешности.

Этот процесс будет корректно работать, только если функция аналитична на области интегрирования; на концах области допустимы устранимые особенности, если предел функции (и всех её производных) существует при приближении к этим концам. Таким образом, romberg можно использовать для интегрирования функции вида f(x)=sin(x)/x на интервале, начинающемся или заканчивающемся в нуле.

Обратите внимание: romberg либо не сойдётся, либо вернёт неверные оценки погрешности, если функция или любая из её производных имеет особенность в какой-либо точке области интегрирования (кроме упомянутого выше случая). Поэтому при интегрировании, например, функции 1/(1-x**2) необходимо избегать точек, в которых производная имеет особенность.

НЕСОБСТВЕННЫЕ ИНТЕГРАЛЫ

Интегрирование методом Ромберга также удобно для вычисления интегралов функций на полунеограниченных интервалах или функций с особенностями. Приём состоит в замене переменной, устраняющей особенность и перемещающей её к одному из концов области интегрирования. Пакет math::calculus предоставляет несколько процедур romberg для распространённых случаев.

romberg_infinity

Вычисляет интеграл функции на полунеограниченном интервале; бесконечным может быть a или b. a и b должны иметь одинаковый знак; если нужно интегрировать через ось, например от отрицательного значения до положительной бесконечности, используйте romberg для интегрирования от отрицательного значения до малого положительного значения, а затем romberg_infinity — для интегрирования от положительного значения до положительной бесконечности. Процедура romberg_infinity выполняет замену переменной u=1/x, поэтому интеграл от a до b функции f(x) вычисляется как интеграл от 1/a до 1/b функции f(1/u)/u**2.
romberg_powerLawLower и romberg_powerLawUpper

Вычисляют интеграл функции, имеющей интегрируемую степенную особенность на нижней или верхней границе области интегрирования (либо производная которой имеет там степенную особенность). Эти процедуры принимают первым параметром gamma, задающий степенной показатель. Предполагается, что функция или её первая производная неограниченно возрастает как (x-a)**(-gamma) или (b-x)**(-gamma). Значение gamma должно быть больше нуля и меньше 1.

Эти процедуры полезны не только для интегрирования функций, стремящихся к бесконечности на одном из концов области, но и функций, производные которых не существуют на её конце. Например, интегрирование f(x)=pow(x,0.25) на области, одним из концов которой является начало координат, приведёт к тому, что процедура romberg значительно занизит погрешность интеграла. Эту проблему можно решить, заметив, что первая производная f(x), f'(x)=x**(-3/4)/4, стремится к бесконечности в начале координат. Интегрирование с помощью romberg_powerLawLower и значением gamma, равным 0.75, обеспечивает гораздо более устойчивую сходимость.

Эти процедуры выполняют замену переменной u=(x-a)**(1-gamma) (romberg_powerLawLower) или u=(b-x)**(1-gamma) (romberg_powerLawUpper).

Итак, значение gamma выбирается следующим образом:
- Если f(x) ~ x**(-a) (0 < a < 1), используйте gamma = a
- Если f'(x) ~ x**(-b) (0 < b < 1), используйте gamma = b
romberg_sqrtSingLower и romberg_sqrtSingUpper

Для распространённого случая gamma=0.5 эти процедуры работают так же, как romberg_powerLawLower и romberg_powerLawUpper: они вычисляют интеграл функции с обратной корневой особенностью на одном из концов интервала. Их реализация проще и использует квадратные корни, а не произвольные степени.
romberg_expLower и romberg_expUpper

Эти процедуры предназначены для интегрирования функции, экспоненциально возрастающей или убывающей на полунеограниченном интервале. romberg_expLower обрабатывает экспоненциально возрастающие функции и допускает сколь угодно большое отрицательное значение нижнего предела интегрирования. romberg_expUpper обрабатывает экспоненциально убывающие функции и допускает сколь угодно большое положительное значение верхнего предела интегрирования. Процедуры выполняют замену переменной u=exp(-x) и u=exp(x) соответственно.

ДРУГИЕ ЗАМЕНЫ ПЕРЕМЕННОЙ

Если вам требуется вычислить несобственный интеграл, не перечисленный здесь, замену переменной можно записать всего в несколько строк Tcl. Поскольку соответствующий код Tcl несколько необычен, здесь приведён подробный пример.

Предположим, что требуется проинтегрировать функцию f(x)=exp(x)/sqrt(1-x*x) (не очень естественная функция, но подходящий пример) на интервале (-1,1). Знаменатель обращается в ноль на обоих концах интервала. Требуется выполнить замену переменной x на u так, чтобы dx/sqrt(1-x**2) преобразовалось в du. Выбрав x=sin(u), получаем dx=cos(u)*du и sqrt(1-x**2)=cos(u). Интеграл от a до b функции f(x) равен интегралу от asin(a) до asin(b) функции f(sin(u))*cos(u).

Можно создать функцию g, принимающую произвольную функцию f и параметр u и вычисляющую новый подынтегральный член.

proc g { f u } {
    set x [expr { sin($u) }]
    set cmd $f; lappend cmd $x; set y [eval $cmd]
    return [expr { $y / cos($u) }]
}

Теперь интегрирование f от a до b эквивалентно интегрированию g от asin(a) до asin(b). Вычислить f в контексте вызывающей процедуры непросто; для этого предназначена следующая процедура.

proc romberg_sine { f a b args } {
    set f [lreplace $f 0 0 [uplevel 1 [list namespace which [lindex $f 0]]]]
    set f [list g $f]
    return [eval [linsert $args 0 romberg $f [expr { asin($a) }] [expr { asin($b) }]]]
}

Эта процедура romberg_sine позволяет интегрировать любую функцию, в знаменателе которой присутствует sqrt(1-x*x). Наша тестовая функция — f(x)=exp(x)/sqrt(1-x*x):

proc f { x } {
    expr { exp($x) / sqrt( 1. - $x*$x ) }
}

Чтобы вычислить её интеграл, достаточно применить romberg_sine так же, как любую другую процедуру romberg:

foreach { value error } [romberg_sine f -1.0 1.0] break
puts [format "integral is %.6g +/- %.6g" $value $error]

integral is 3.97746 +/- 2.3557e-010

Ошибки, идеи, отзывы

В этом документе и описываемом в нём пакете наверняка есть ошибки и другие проблемы. Пожалуйста, сообщайте о них в категории math :: calculus на сайте системы отслеживания Tcllib. Также сообщайте о любых идеях по улучшению пакета и/или документации.

Предлагая изменения кода, пожалуйста, прикладывайте унифицированные различия, то есть вывод команды diff -u.

Также обратите внимание, что предпочтительнее прикладывать вложения, а не вставлять исправления непосредственно в текст. Чтобы добавить вложение, сразу после создания заявки откройте форму Edit и нажмите самую левую кнопку на дополнительной панели навигации.

СМ. ТАКЖЕ

math::calculus, math::interpolate

КАТЕГОРИЯ

Математика

АВТОРСКИЕ ПРАВА

Licensed under the BSD license
https://core.tcl-lang.org/tcllib/doc/trunk/embedded/md/tcllib/files/modules/math/romberg.md