23.1 Функции одной переменной ¶

Octave поддерживает пять различных адаптивных алгоритмов квадратур для вычисления интеграла функции f на отрезке от a до b. Это

quad

Численное интегрирование, основанное на квадратуре Гаусса.

quadv

Численное интегрирование с помощью адаптивного векторизованного правила Симпсона.

quadl

Численное интегрирование с помощью адаптивного правила Лобатто.

quadgk

Численное интегрирование с помощью адаптивного правила Гаусса-Конрода.

quadcc

Численное интегрирование с помощью адаптивных правил Кленшоу-Картиса.

Кроме того, также предоставлены следующие функции:

integral

Функция-обёртка совместимости, которая будет выбирать между quadv и quadgk в зависимости от подынтегральной функции и выбранных параметров.

trapz, cumtrapz

Численное интегрирование данных с помощью метода трапеций.

Лучший алгоритм квадратур зависит от подынтегральной функции. Если у вас есть эмпирические данные, а не функция, выбор — trapz или cumtrapz. Если вы не уверены в характеристиках подынтегральной функции, quadcc будет наиболее надёжным, так как он может обрабатывать разрывы, особенности, осциллирующие функции и бесконечные промежутки. Когда подынтегральная функция гладкая, quadgk может быть самым быстрым из алгоритмов.

	Функция	Характеристики
	quad	Низкая точность при негладких подынтегральных функциях
	quadv	Средняя точность при гладких подынтегральных функциях
	quadl	Средняя точность при гладких подынтегральных функциях. Медленнее, чем quadgk.
	quadgk	Средняя точность (1e-6 – 1e-9) при гладких подынтегральных функциях.
		Обрабатывает осциллирующие функции и бесконечные пределы
	quadcc	Низкая до высокой точности при негладких/гладких подынтегральных функциях
		Обрабатывает осциллирующие функции, особенности и бесконечные пределы

Вот пример использования quad для интегрирования функции

f(x) = x * sin (1/x) * sqrt (abs (1 - x))

от x = 0 до x = 3.

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

Первый шаг — определение функции:

function y = f (x)
  y = x .* sin (1./x) .* sqrt (abs (1 - x));
endfunction

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

Второй шаг — вызов quad с пределами интегрирования:

[q, ier, nfun, err] = quad ("f", 0, 3)
     ⇒ 1.9819
     ⇒ 1
     ⇒ 5061
     ⇒ 1.1522e-07

Хотя quad возвращает ненулевое значение для ier, результат достаточно точен (чтобы увидеть почему, посмотрите, что произойдёт с результатом, если вы сместите нижний предел до 0.1, затем 0.01, затем 0.001 и т. д.).

Функция "f" может быть именем функции или функцией-обработчиком. Эти варианты делают интегрирование довольно простым, не прибегая к полному определению функции в файле m-формата. Например:

# Verify gamma function = (n-1)! for n = 4
f = @(x) x.^3 .* exp (-x);
quadcc (f, 0, Inf)
     ⇒ 6.0000

: q = quad (f, a, b) ¶

: q = quad (f, a, b, tol) ¶

: q = quad (f, a, b, tol, sing) ¶

: [q, ier, nfev, err] = quad (…) ¶

Численно вычислить интеграл от f от a до b с использованием Fortran-подпрограмм из QUADPACK.

f — это функция-обработчик, встроенная функция или строка, содержащая имя функции для вычисления. Функция должна иметь вид y = f (x), где y и x — скаляры.

a и b — нижний и верхний пределы интегрирования. Один или оба могут быть бесконечными.

Необязательный аргумент tol — вектор, который определяет желаемую точность результата. Первый элемент вектора — желаемая абсолютная погрешность, а второй — желаемая относительная погрешность. Для выбора только относительного теста установите абсолютную погрешность в ноль. Для выбора только абсолютного теста установите относительную погрешность в ноль. Обе погрешности по умолчанию равны sqrt (eps) или приблизительно 1.5e-8.

Необязательный аргумент sing — вектор значений, в которых подынтегральная функция известна как особенная.

Результат интегрирования возвращается в q.

ier содержит целочисленный код ошибки (0 указывает на успешное интегрирование).

nfev указывает количество вычислений функции.

err содержит оценку погрешности решения.

Функция quad_options может установить другие необязательные параметры для quad.

Примечание: поскольку quad написана на Fortran, она не может вызываться рекурсивно. Это предотвращает её использование в интегрировании по нескольким переменным подпрограммами dblquad и triplequad.

См. также: quad_options, quadv, quadl, quadgk, quadcc, trapz, dblquad, triplequad.

: quad_options () ¶

: val = quad_options (opt) ¶

: quad_options (opt, val) ¶

Запрос или установка параметров для функции quad.

При вызове без аргументов отображаются имена всех доступных параметров и их текущие значения.

При вызове с одним аргументом возвращает значение параметра opt.

При вызове с двумя аргументами quad_options устанавливает параметр opt в значение val.

Параметры включают

"absolute tolerance"

Абсолютная погрешность; может быть равна нулю для чистого относительного теста.

"relative tolerance"

Неотрицательная относительная погрешность. Если абсолютная погрешность равна нулю, относительная погрешность должна быть больше или равна max (50*eps, 0.5e-28).

"single precision absolute tolerance"

Абсолютная погрешность для одинарной точности; может быть равна нулю для чистого относительного теста.

"single precision relative tolerance"

Неотрицательная относительная погрешность для одинарной точности. Если абсолютная погрешность равна нулю, относительная погрешность должна быть больше или равна max (50*eps, 0.5e-28).

: q = quadv (f, a, b) ¶

: q = quadv (f, a, b, tol) ¶

: q = quadv (f, a, b, tol, trace) ¶

: q = quadv (f, a, b, tol, trace, p1, p2, …) ¶

: [q, nfev] = quadv (…) ¶

Числовое вычисление интеграла от f от a до b с использованием адаптивного правила Симпсона.

f — это функция-обработчик, встроенная функция или строка, содержащая имя функции для вычисления. quadv является векторизованной версией quad и функция, определённая f, должна принимать скаляр или вектор в качестве входных данных и возвращать скаляр, вектор или массив в качестве выходных данных.

a и b — нижний и верхний пределы интегрирования. Оба предела должны быть конечными.

Необязательный аргумент tol определяет абсолютную толерантность, используемую для остановки процедуры адаптации. Значение по умолчанию равно 1e-6.

Алгоритм, используемый quadv, включает рекурсивное деление интервала интегрирования и применение правила Симпсона на каждом подынтервале. Если trace имеет значение true, то после вычисления каждого из этих частных интегралов отображаются: (1) общее количество вычислений функции, (2) левая граница подынтервала, (3) длина подынтервала, (4) приближение интеграла по подынтервалу.

Дополнительные аргументы p1 и т. д. передаются непосредственно функции f. Для использования значений по умолчанию для tol и trace можно передать пустые матрицы ([]).

Результат интегрирования возвращается в q.

Необязательный выход nfev указывает общее количество вычислений функции.

Примечание: quadv написана на языке сценариев Octave и может использоваться рекурсивно в dblquad и triplequad, в отличие от функции quad.

См. также: quad, quadl, quadgk, quadcc, trapz, dblquad, triplequad, integral, integral2, integral3.

: q = quadl (f, a, b) ¶

: q = quadl (f, a, b, tol) ¶

: q = quadl (f, a, b, tol, trace) ¶

: q = quadl (f, a, b, tol, trace, p1, p2, …) ¶

: [q, nfev] = quadl (…) ¶

Числовое вычисление интеграла от f от a до b с использованием адаптивного правила Лобатто.

f — функция-обработчик, встроенная функция или строка, содержащая имя функции для вычисления. Функция f должна быть векторизованной и возвращать вектор выходных значений при получении вектора входных значений.

a и b — нижний и верхний пределы интегрирования. Оба предела должны быть конечными.

Необязательный аргумент tol определяет абсолютную толерантность, с которой выполняется интегрирование. Значение по умолчанию равно 1e-6.

Алгоритм, используемый quadl, включает рекурсивное деление интервала интегрирования. Если trace определено, то для каждого подынтервала отображаются: (1) общее количество вычислений функции, (2) левая граница подынтервала, (3) длина подынтервала, (4) приближение интеграла по подынтервалу.

Дополнительные аргументы p1 и т. д. передаются непосредственно функции f. Для использования значений по умолчанию для tol и trace можно передать пустые матрицы ([]).

Результат интегрирования возвращается в q.

Необязательный выход nfev указывает общее количество вычислений функции.

Ссылка: W. Gander и W. Gautschi, Adaptive Quadrature - Revisited, BIT Vol. 40, No. 1, March 2000, pp. 84–101. https://www.inf.ethz.ch/personal/gander/

См. также: quad, quadv, quadgk, quadcc, trapz, dblquad, triplequad, integral, integral2, integral3.

: q = quadgk (f, a, b) ¶

: q = quadgk (f, a, b, abstol) ¶

: q = quadgk (f, a, b, abstol, trace) ¶

: q = quadgk (f, a, b, "prop", val, …) ¶

: [q, err] = quadgk (…) ¶

Численно вычисляет интеграл функции f от a до b с помощью адаптивного квадратурного метода Гаусса-Кронерода.

f — это дескриптор функции, встроенная функция или строка, содержащая имя функции для вычисления. Функция f должна быть векторизованной и возвращать вектор выходных значений при вводе вектора входных значений (см. свойство "ArrayValued" для исключения из этого правила).

a и b — нижний и верхний пределы интегрирования. Один или оба предела могут быть бесконечными или содержать слабые особые точки на концах. Будет использовано преобразование переменных для обработки бесконечных интервалов и ослабления особенностей. Например:

quadgk (@(x) 1 ./ (sqrt (x) .* (x + 1)), 0, Inf)

Обратите внимание, что формулировка подынтегральной функции использует оператор поэлементного умножения ./ , и все пользовательские функции для quadgk должны делать то же самое.

Необязательный аргумент abstol определяет абсолютную погрешность, используемую для остановки процедуры интегрирования. Значение по умолчанию — 1e-10 (1e-5 для типа single).

Алгоритм, используемый функцией quadgk, включает разделение интервала интегрирования и вычисление каждого подынтервала. Если trace имеет значение true, то после вычисления каждого из этих частичных интегралов отображаются: (1) количество подынтервалов на данном шаге, (2) текущая оценка погрешности err, (3) текущая оценка интеграла q.

Поведение алгоритма можно настроить, передавая аргументы в quadgk в виде пар "prop", val . Допустимые свойства:

AbsTol

Определяет абсолютную погрешность для квадратурного метода. Значение по умолчанию — 1e-10 (1e-5 для типа single).

RelTol

Определяет относительную погрешность для квадратурного метода. Значение по умолчанию — 1e-6 (1e-4 для типа single).

ArrayValued

При установке в значение true, функция f возвращает массив для скалярного входного значения. По умолчанию значение false требует, чтобы функция f возвращала выходные данные того же размера, что и входные данные. Например,

quadgk (@(x) x .^ (1:5), 0, 2, "ArrayValued", 1)

вычислит интеграл [x.^1, x.^2, x.^3, x.^4, x.^5] в одном вызове функции, вместо необходимости многократно определять анонимную функцию и использовать обычный вызов quadgk.

WayPoints

Указывает точки, которые станут границами подынтервалов в алгоритме, что может привести к значительно улучшенной оценке погрешности интеграла, более быстрому вычислению или обоим. Это может быть полезно для указания дополнительных подынтервалов вокруг области, где подынтегральная функция быстро изменяется, или для обозначения мест разрыва первой производной функции. Например, функция signum имеет разрыв в x == 0, и путем указания контрольной точки

quadgk (@(x) sign (x), -0.5, 1, "Waypoints", [0])

предел погрешности уменьшается с 4e-7 до 1e-13.

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

Если какие-либо из контрольных точек являются комплексными, выполняется контурное интегрирование, как описано ниже.

MaxIntervalCount

quadgk изначально делит интервал для выполнения квадратурного метода на 10 интервалов или, если заданы контрольные точки, в каждой контрольной точке. Подынтервалы, которые имеют неприемлемую погрешность, делятся и перевычисляются. Если количество подынтервалов превышает 650 в какой-либо момент, то сигнализируется о плохой сходимости, и возвращается текущая оценка интеграла. Свойство "MaxIntervalCount" может быть использовано для изменения количества подынтервалов, которое может существовать перед выходом.

Trace

Если логически истинно quadgk, выводится информация о сходимости квадратурного метода на каждом шаге.

Если какие-либо из a, b или waypoints являются комплексными, квадратурный метод рассматривается как контурный интеграл вдоль кусочно-линейной траектории, определенной [a, waypoints(1), waypoints(2), …, b]. В этом случае предполагается, что интеграл не имеет особых точек на границе. Например,

quadgk (@(z) log (z), 1+1i, 1+1i, "WayPoints",
        [-1+1i, -1-1i, +1-1i])

вычисляет интеграл log (z) вдоль квадрата, определяемого [1+1i, -1+1i, -1-1i, +1-1i].

Результат интегрирования возвращается в q.

err — приблизительная оценка погрешности интеграла abs (q - I), где I — точное значение интеграла. Если адаптивное интегрирование не сошлось, значение err будет больше, чем заданная точность. Если требуется только один выход, то будет выведено предупреждение, когда заданная точность не достигнута. Если запрашиваются оба выхода, то предупреждение не выдается, и программист отвечает за проверку и определение того, являются ли результаты удовлетворительными.

Ссылка: Л.Ф. Шэмпин, ""Векторизованная адаптивная квадратура в MATLAB"", Журнал вычислительной и прикладной математики, стр. 131–140, том 211, выпуск 2, февраль 2008.

См. также: quad, quadv, quadl, quadcc, trapz, dblquad, triplequad, integral, integral2, integral3.

: q = quadcc (f, a, b) ¶

: q = quadcc (f, a, b, tol) ¶

: q = quadcc (f, a, b, tol, sing) ¶

: [q, err, nr_points] = quadcc (…) ¶

Численно вычислить интеграл функции f от a до b с использованием двукратно адаптивной квадратуры Кленшоу-Кёртиса.

f — это функция-обработчик, встроенная функция или строка, содержащая имя функции для вычисления. Функция f должна быть векторизованной и должна возвращать вектор выходных значений, если ей передаётся вектор входных значений. Например,

f = @(x) x .* sin (1./x) .* sqrt (abs (1 - x));

что использует элемент-по-элементную «точечную» форму для всех операторов.

a и b — это нижний и верхний пределы интегрирования. Один или оба предела могут быть бесконечными. quadcc обрабатывает бесконечный предел, заменяя переменную интегрирования на x = tan (pi/2*u).

Необязательный аргумент tol — это вектор из 1 или 2 элементов, который задаёт требуемую точность результата. Первый элемент вектора — желаемая абсолютная погрешность, а второй — желаемая относительная погрешность. Чтобы выбрать только относительный тест, установите абсолютную погрешность в ноль. Чтобы выбрать только абсолютный тест, установите относительную погрешность в ноль. По умолчанию абсолютная погрешность равна 1e-10 (1e-5 для single), а относительная погрешность — 1e-6 (1e-4 для single).

Необязательный аргумент sing содержит список точек, в которых подынтегральная функция имеет известные особенности или разрывы в любой из её производных внутри интервала интегрирования. Для приведённого выше примера, имеющего разрыв в x=1, вызов quadcc будет следующим

int = quadcc (f, a, b, [], [ 1 ]);

Результат интегрирования возвращается в q.

err — это оценка абсолютной погрешности интегрирования.

nr_points — это количество точек, в которых была вычислена подынтегральная функция.

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

quadcc способен обрабатывать нечисловые значения подынтегральной функции, такие как NaN или Inf. Если интеграл расходится и quadcc обнаруживает это, выдаётся предупреждение, и возвращается Inf или -Inf.

Примечание: quadcc — это универсальный алгоритм квадратуры и, как таковой, может быть менее эффективным для гладкой или иначе хорошо себя ведущей подынтегральной функции, чем другие методы, такие как quadgk.

Алгоритм использует правила квадратуры Кленшоу-Кёртиса возрастающей степени в каждом интервале и делит интервал пополам, если функция не кажется гладкой или достигнута степень максимальной степени. Оценка погрешности вычисляется из L2-нормы разности между двумя последовательными интерполяциями подынтегральной функции по узлам соответствующих правил квадратуры.

Ссылка: P. Gonnet, Увеличение надёжности адаптивной квадратуры с использованием явных интерполянтов, ACM Transactions on Mathematical Software, Vol. 37, Issue 3, Article No. 3, 2010.

См. также: quad, quadv, quadl, quadgk, trapz, dblquad, triplequad.

: q = integral (f, a, b) ¶

: q = integral (f, a, b, prop, val, …) ¶

: [q, err] = integral (…) ¶

Численно вычислить интеграл функции f от a до b с использованием адаптивной квадратуры.

integral — это оболочка для quadcc (общие вещественные скалярные подынтегральные функции и пределы) и quadgk (интегралы со специфицированными путями интегрирования и подынтегральными функциями с массивами значений), предназначенная для обеспечения совместимости с MATLAB. Более тонкий контроль численного интегрирования может быть достигнут путём вызова различных функций квадратуры напрямую.

f — это функция-обработчик, встроенная функция или строка, содержащая имя функции для вычисления. Функция f должна быть векторизованной и возвращать вектор выходных значений, если ей передаётся вектор входных значений.

a и b — это нижний и верхний пределы интегрирования. Один или оба предела могут быть бесконечными или содержать слабые особенности на концах. Если один или оба предела являются комплексными, integral будет выполнять интеграл по прямой. В качестве альтернативы, путь сложной области можно указать с помощью опции "Waypoints" (см. ниже).

Дополнительные необязательные параметры могут быть заданы с помощью пар "property", value. Допустимые свойства:

Waypoints

Определяет точки, используемые для определения подынтервалов алгоритма квадратуры, или если a, b или waypoints являются комплексными, то квадратура вычисляется как контурный интеграл вдоль кусково-непрерывного пути. Подробнее см. quadgk.

ArrayValued

integral ожидает, что f вернёт скалярное значение, если не задано arrayvalued = true. Этот параметр заставит integral выполнить интегрирование по всему массиву и вернуть q с теми же размерностями, что и возвращает f. Подробнее см. quadgk.

AbsTol

Определите абсолютную погрешность для квадратуры. По умолчанию абсолютная погрешность равна 1e-10 (1e-5 для single).

RelTol

Определите относительную погрешность для квадратуры. По умолчанию относительная погрешность равна 1e-6 (1e-4 для single).

Необязательный выходной параметр err содержит оценку абсолютной погрешности, используемой вызываемым интегрирователем.

Адаптивная квадратура используется для минимизации оценки погрешности до тех пор, пока не выполняется следующее условие:

error <= max (AbsTol, RelTol*|q|).

Известные несовместимости с MATLAB:

1. Если толерантности не указаны, а любые пределы интегрирования или контрольные точки имеют тип single, то функции интегрирования Octave автоматически уменьшают значения по умолчанию абсолютной и относительной погрешностей, как указано выше. Если требуются более строгие толерантности, они должны быть указаны. MATLAB оставляет более строгие толерантности, соответствующие входам типа double неизменными независимо от типа пределов интегрирования.

См. также: integral2, integral3, quad, quadgk, quadv, quadl, quadcc, trapz, dblquad, triplequad.

Иногда у нас нет функции, а есть только сырые точки (x, y), по которым необходимо выполнить интегрирование. Это может произойти при сборе данных в эксперименте. Функция trapz может интегрировать эти значения, как показано в следующем примере, где данные "data" были собраны по функции косинуса на интервале [0, pi/2).

x = 0:0.1:pi/2;  # Uniformly spaced points
y = cos (x);
trapz (x, y)
     ⇒ 0.99666

Ответ довольно близок к точному значению 1. Обычная квадратура чувствительна к характеристикам подынтегральной функции. Эмпирическое интегрирование зависит не только от подынтегральной функции, но и от выбранных для представления функции точек. Повторение приведённого выше примера с функцией синуса на интервале [0, pi/2) даёт гораздо худшие результаты.

x = 0:0.1:pi/2;  # Uniformly spaced points
y = sin (x);
trapz (x, y)
     ⇒ 0.92849

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

x = linspace (0, pi/2, 16);  # Uniformly spaced, but including endpoint
y = sin (x);
trapz (x, y)
     ⇒ 0.99909

В общем случае может не быть способа узнать наилучшее распределение точек заранее. Или точки могут поступать из эксперимента, где нет возможности выбрать наилучшее распределение. В любом случае необходимо учитывать этот вопрос при использовании trapz.

: q = trapz (y) ¶

: q = trapz (x, y) ¶

: q = trapz (…, dim) ¶

Численно вычислить интеграл точек y с использованием метода трапеций.

trapz (y) вычисляет интеграл y по первому неединичному измерению. Если аргумент x опущен, предполагается равномерно распределённый вектор x с шагом 1. trapz (x, y) вычисляет интеграл относительно шага в x и значений в y. Это полезно, если точки в y были взяты неравномерно.

Если задан необязательный аргумент dim, то работа производится вдоль этого измерения.

Примечание: Если x не указан, будет использован шаг 1. Для масштабирования интеграла до правильного значения необходимо умножить на фактический шаг (deltaX). Например, интеграл от x^3 в пределах от [0, 1] равен x^4/4 или 0,25. Следующий код использует trapz для вычисления интеграла тремя различными способами.

x = 0:0.1:1;
y = x.^3;
## No scaling
q = trapz (y)
  ⇒ q = 2.5250
## Approximation to integral by scaling
q * 0.1
  ⇒ 0.25250
## Same result by specifying x
trapz (x, y)
  ⇒ 0.25250

См. также: cumtrapz.

: q = cumtrapz (y) ¶

: q = cumtrapz (x, y) ¶

: q = cumtrapz (…, dim) ¶

Кумулятивное численное интегрирование точек y с использованием метода трапеций.

cumtrapz (y) вычисляет кумулятивный интеграл y по первому неединичному измерению. В то время как trapz сообщает только сумму интеграла, cumtrapz сообщает текущее значение частичной суммы в каждой точке y.

Если аргумент x опущен, предполагается равномерно распределённый вектор x с шагом 1. cumtrapz (x, y) вычисляет интеграл относительно шага в x и значений в y. Это полезно, если точки в y были взяты неравномерно.

Если задан необязательный аргумент dim, то работа производится вдоль этого измерения.

Примечание: Если x не указан, будет использован шаг 1. Для масштабирования интеграла до правильного значения необходимо умножить на фактический шаг (deltaX).

См. также: trapz, cumsum.