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, nfun, 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 указывает на успешное интегрирование).

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

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, nfun] = quadv (…) ¶

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

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

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

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

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

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

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

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

Примечание: 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, nfun] = quadl (…) ¶

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

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

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

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

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

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

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

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

Ссылка: W. Gander и W. Gautschi, Адаптивное квадратурное интегрирование - пересмотр, BIT Vol. 40, No. 1, март 2000, стр. 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 должна быть векторизованной и возвращать вектор выходных значений при передаче вектора входных значений.

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

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

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

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

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

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

AbsTol

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

RelTol

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

MaxIntervalCount

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

WayPoints

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

quadgk (@(x) abs (1 - x.^2), 0, 2, "Waypoints", 1)

указывают точку разрыва в подынтегральной функции в x = 1.

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 — точное значение интеграла.

Ссылка: Л.Ф. Шампин, "Векторизованная адаптивная квадратура в 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 для одинарной точности), а значение относительной погрешности по умолчанию — 1e-6 (1e-4 для одинарной точности).

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

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

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

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

nr_points — количество точек, в которых вычислялось подынтегральное выражение.

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

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

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

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

Ссылка: П. Гоннет, Повышение надежности адаптивной квадратуры с использованием явных интерполянтов, ACM Transactions on Mathematical Software, том 37, выпуск 3, статья № 3, 2010 г.

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

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

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

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

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

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

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

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

Waypoints

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

ArrayValued

integral ожидает, что f вернет скалярное значение, если не задано arrayvalued = true. Этот параметр заставит integral выполнить интегрирование по всему массиву и вернуть q с теми же размерностями, что и возвращает f. Дополнительную информацию см. в quadv.

AbsTol

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

RelTol

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

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

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

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

1. Если погрешности не указаны, а любые пределы интегрирования или точки пути являются типа single, то функции интегрирования Octave автоматически уменьшают значения по умолчанию для абсолютной и относительной погрешностей, как указано выше. Если требуется более высокая точность, необходимо указать значения погрешностей. MATLAB оставляет более строгие значения погрешности, соответствующие входным данным типа double независимо от типа пределов интегрирования.
2. Вследствие использования quadcc, quadgk, и quadv, определённые сочетания опций не поддерживаются. В настоящее время "ArrayValued" не может быть использовано совместно с "RelTol" или "Waypoints".

См. также: 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 не указан, используется единичный шаг. Чтобы правильно масштабировать интеграл, необходимо умножить на фактический шаг (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 не указан, используется единичный шаг. Для корректного масштабирования интеграла необходимо умножить на фактический шаг (deltaX).

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