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

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

quad

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

quadv

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

quadl

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

quadgk

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

quadcc

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

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

integral

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

trapz, cumtrapz

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

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

Вот пример использования 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 integral (x^3) = x^4/4
f = inline ("x.^3");
quadgk (f, 0, 1)
     ⇒ 0.25000

# 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 — true, то после вычисления каждого из этих частичных интегралов отображаются: (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 том 40, № 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)

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

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

Используемый алгоритм quadgk включает деление интервала интегрирования и вычисление каждого подынтервала. Если trace истинно, то после вычисления каждого из этих частичных интегралов отображаются: (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 является комплексным, то квадратура обрабатывается как контурный интеграл по кусочно-непрерывной траектории, определенной выше. В этом случае предполагается, что интеграл не имеет особенностей на границе. Например,

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

Ссылка: L.F. Shampine, «Векторизованное адаптивное квадратурное интегрирование в 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-нормы разности между двумя последовательными интерполяциями подынтегральной функции по узлам соответствующих правил квадратуры.

Примечание реализации: для версий Octave ≤ 4.2 quadcc принимал единственный аргумент погрешности, который определял относительную погрешность. Для версий 4.4 и 5 quadcc выведет предупреждение при вызове с единственным аргументом погрешности, указывая, что значение этого входа изменилось с относительной погрешности на абсолютную погрешность. Идентификатор предупреждения для этого сообщения — "Octave:quadcc:RelTol-conversion". Предупреждение можно отключить с помощью warning ("off", "Octave:quadcc:RelTol-conversion").

Ссылка: P. Gonnet, Повышение надежности адаптивного квадратурного интегрирования с использованием явных интерполянтов, 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. Если погрешности не указаны, а любые пределы интегрирования или waypoints имеют тип single, то функции интегрирования Octave автоматически уменьшают значения по умолчанию абсолютной и относительной погрешностей, как указано выше. Если требуются более узкие допуски, они должны быть указаны. MATLAB оставляет более строгие допуски, соответствующие double входным данным, независимо от типа пределов интегрирования.
2. Вследствие использования quadcc, quadgk, и quadv, некоторые комбинации параметров не поддерживаются. В настоящее время "ArrayValued" нельзя комбинировать с "RelTol" или "Waypoints".

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

Иногда у вас нет функции, а есть только сырые точки (x, y), по которым нужно выполнить интегрирование. Это может произойти при сборе данных в эксперименте. Функция trapz может интегрировать эти значения, как показано в следующем примере, где данные были собраны по функции косинуса в диапазоне [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.