Открыть сервис

Адаптивная квадратура Гаусса-Кронрода

Адаптивная квадратура Гаусса-Кронрода — это численный метод приближённого вычисления определённых интегралов, сочетающий квадратурные формулы Гаусса и Кронрода с автоматическим выбором шага интегрирования (адаптацией) для достижения заданной точности. Метод широко применяется в вычислительной математике, физике и инженерных расчётах благодаря высокой точности и возможности контроля погрешности.

История

Квадратура Гаусса была разработана Карлом Фридрихом Гауссом в начале XIX века. Она основана на выборе узлов интегрирования как корней ортогональных многочленов (например, многочленов Лежандра) и обеспечивает максимальную алгебраическую точность для заданного числа узлов. Однако недостатком метода является отсутствие простой оценки погрешности, так как для её вычисления требуется дополнительное увеличение числа узлов, что приводит к пересчёту всех значений подынтегральной функции.

В 1964 году Александр Кронрод предложил модификацию квадратуры Гаусса, добавив к узлам Гаусса дополнительные узлы, расположенные между ними. Это позволило получить квадратурную формулу с более высокой точностью и одновременно оценить погрешность интегрирования. В 1970-х годах идея Кронрода была объединена с адаптивными алгоритмами, что привело к созданию эффективных методов численного интегрирования, реализованных в библиотеках численного анализа (например, QUADPACK).

Основные принципы

Квадратура Гаусса-Лежандра

Квадратура Гаусса-Лежандра для интеграла $\int_{-1}^{1} f(x) \, dx$ использует $n$ узлов $x_i$ (корни многочлена Лежандра $P_n(x)$) и веса $w_i$, определяемые по формуле: \[ \int_{-1}^{1} f(x) \, dx \approx \sum_{i=1}^{n} w_i f(x_i). \] Эта формула точна для многочленов степени не выше $2n-1$. Для интеграла на произвольном отрезке $[a, b]$ выполняется линейное преобразование.

Квадратура Кронрода

Кронрод предложил добавить к $n$ узлам Гаусса ещё $n+1$ узел, расположенный между ними, включая концы отрезка. Общее число узлов становится $2n+1$. Веса для новых узлов вычисляются так, чтобы формула была точна для многочленов степени не выше $3n+1$ (при $n$ нечётном) или $3n$ (при $n$ чётном). Например, для $n=7$ (7 узлов Гаусса) используется 15 узлов Кронрода (7+8). Погрешность оценивается как разность между результатами квадратур Гаусса и Кронрода: \[ E \approx \left| I_{\text{Кронрод}} - I_{\text{Гаусс}} \right|. \]

Адаптивный алгоритм

Адаптивная квадратура Гаусса-Кронрода разбивает отрезок интегрирования на подотрезки, на каждом из которых применяется пара квадратур (Гаусса и Кронрода). Если оценка погрешности на текущем подотрезке превышает заданный допуск (абсолютный или относительный), отрезок делится пополам, и процедура повторяется рекурсивно. Алгоритм продолжается до тех пор, пока на всех подотрезках не будет достигнута требуемая точность. Это позволяет эффективно обрабатывать функции с резкими изменениями, осцилляциями или особенностями.

Реализация

Алгоритм (псевдокод)

`` function adaptive_gauss_kronrod(f, a, b, tol, max_depth): I_gauss, I_kronrod = compute_gauss_kronrod(f, a, b) err = abs(I_kronrod - I_gauss) if err < tol or depth >= max_depth: return I_kronrod else: mid = (a + b) / 2 return adaptive_gauss_kronrod(f, a, mid, tol/2, depth+1) + adaptive_gauss_kronrod(f, mid, b, tol/2, depth+1) ``

Популярные реализации

  • QUADPACK (Fortran) — библиотека, разработанная в 1970-х годах, включает подпрограмму dqag (адаптивная квадратура Гаусса-Кронрода). Она используется в GNU Scientific Library (GSL) и SciPy.
  • SciPy (Python) — функция scipy.integrate.quad основана на QUADPACK и использует адаптивную квадратуру Гаусса-Кронрода с 15, 21, 31, 41, 51 или 61 узлом.
  • MATLAB — функция integral (с 2012 года) реализует адаптивную квадратуру Гаусса-Кронрода с автоматическим выбором числа узлов.
  • GSL (C/C++) — функция gsl_integration_qag предоставляет адаптивный алгоритм с несколькими правилами Кронрода.

Примеры применения

Пример 1: Гладкая функция

Вычисление интеграла $\int_{0}^{1} e^{-x^2} \, dx$:

  • Точное значение: $\frac{\sqrt{\pi}}{2} \text{erf}(1) \approx 0.746824132812427$.
  • Адаптивная квадратура Гаусса-Кронрода (15 узлов) даёт результат с точностью $10^{-12}$ за 1-2 итерации, так как функция гладкая.

Пример 2: Функция с особенностью

Вычисление $\int_{0}^{1} \frac{1}{\sqrt{x}} \, dx = 2$:

  • Подынтегральная функция имеет особенность в нуле. Адаптивный алгоритм автоматически сгущает узлы вблизи $x=0$, достигая точности $10^{-10}$ за 10-15 делений.

Пример 3: Осциллирующая функция

Вычисление $\int_{0}^{10} \sin(100x) \, dx$:

  • Высокая частота осцилляций требует большого числа узлов. Адаптивная квадратура Гаусса-Кронрода эффективно разбивает отрезок на мелкие подотрезки, где функция ведёт себя как синусоида.

Преимущества и недостатки

Преимущества

  • Высокая точность: квадратура Кронрода обеспечивает алгебраическую точность до $3n+1$ степени, что превосходит обычную квадратуру Гаусса.
  • Автоматический контроль погрешности: оценка погрешности по разности двух квадратур не требует дополнительных вычислений.
  • Адаптивность: алгоритм автоматически сгущает узлы в областях с большими градиентами, что экономит вычислительные ресурсы.
  • Робастность: метод хорошо работает для большинства гладких и кусочно-гладких функций.

Недостатки

  • Вычислительная сложность: на каждом шаге требуется вычисление подынтегральной функции во всех узлах (до $2n+1$ точек), что может быть затратно для сложных функций.
  • Чувствительность к разрывам: для функций с разрывами первого рода или сингулярностями алгоритм может потребовать очень мелкого разбиения или не сойтись.
  • Ограничение на размерность: метод применим только к одномерным интегралам; для многомерных интегралов требуются другие подходы (например, кубатурные формулы или методы Монте-Карло).

Сравнение с другими методами

МетодТочностьАдаптивностьОценка погрешностиСложность
Квадратура ГауссаВысокая (для гладких функций)НетНетНизкая
Квадратура КронродаОчень высокаяНетДаСредняя
Адаптивная квадратура Гаусса-КронродаОчень высокаяДаДаВысокая
Метод СимпсонаСредняяДаДаНизкая
Метод Монте-КарлоНизкая (сходимость $1/\sqrt{N}$)НетСтатистическаяВысокая

Применение

  • Физика: вычисление потенциалов, интегралов перекрытия в квантовой механике, интегралов от функций Грина.
  • Инженерия: расчёт моментов инерции, центров масс, электрических полей.
  • Статистика: численное интегрирование плотностей вероятностей, вычисление моментов распределений.
  • Финансовая математика: оценка стоимости опционов (например, метод Блэка-Шоулза требует численного интегрирования).
  • Численное решение дифференциальных уравнений: метод Галёркина и метод конечных элементов часто требуют вычисления интегралов от базисных функций.

Интересные факты

  • Квадратура Кронрода была разработана для улучшения квадратуры Гаусса, но сам Кронрод не использовал адаптивный алгоритм; адаптивность была добавлена позже в библиотеке QUADPACK.
  • В библиотеке QUADPACK реализовано несколько правил Кронрода с разным числом узлов: 15, 21, 31, 41, 51 и 61. Обычно используется 15 узлов (7 Гаусса + 8 Кронрода), так как это даёт хороший баланс между точностью и скоростью.
  • Адаптивная квадратура Гаусса-Кронрода является стандартным методом для одномерного интегрирования в большинстве научных вычислительных пакетов (SciPy, MATLAB, GSL, NAG).
  • Существуют модификации метода для работы с функциями, заданными таблично, или для интегралов с бесконечными пределами (например, преобразование к конечному отрезку).

Источники

  • Кронрод, А. С. (1964). «О приближённом вычислении определённых интегралов». Доклады АН СССР, 154(2), 283–286.
  • Piessens, R., de Doncker-Kapenga, E., Überhuber, C. W., & Kahaner, D. K. (1983). QUADPACK: A Subroutine Package for Automatic Integration. Springer-Verlag.
  • Press, W. H., Teukolsky, S. A., Vetterling, W. T., & Flannery, B. P. (2007). Numerical Recipes: The Art of Scientific Computing (3rd ed.). Cambridge University Press.
  • Документация SciPy: scipy.integrate.quad (https://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.quad.html).

BFOmetr — база данных и аналитика по компаниям России.

На главную BFOmetr →