Ctrl+K

Квадратуры Гаусса

Квадратуры Гаусса#

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

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

В квадратурах Ньютона-Гаусса мы «варьировали» только коэффициенты \(w_i\), тогда как сейчас поиграемся и с положениями узлов \(x_i\).

Рассмотрим интегралы на отрезке \([-1, 1]\):

\[ I=\int_{-1}^1 f(x) d x \approx I_G=\sum_{i=0}^{n-1} w_i f\left(x_i\right) \]

Легко понять, что любая задача может быть сведена к данному виду заменой

\[ \xi=\frac{1}{2}(b+a)+\frac{1}{2}(b-a) x \]
\[ \int_a^b f(\xi) d \xi=\left(\frac{b-a}{2}\right) \int_{-1}^1 f(\xi(x)) d x \]

Используя многочлены Лежандра можно вывести аналитически положения узлов \(x_i\) и весов \(w_i\), которые подойдут для «любой» заранее заданной функции. В данном случае, \(x_i\) будут корнями многочленов Лежандра (подробнее см. квадратуры Гаусса-Кристофеля)

\[\begin{split} \begin{array}{cccc} \hline \begin{array}{c} \text { Число точек } n \\ \text { в квадратуре } \end{array} & \begin{array}{c} \text { Весовые } \\ \text { коэффициенты } \end{array} & \begin{array}{c} \text { Значения } \\ \text { аргумента } \end{array} & \begin{array}{c} \text { Погрешность } \\ \text { аппроксимации } \end{array} \\ \hline 2 & w_0=1,000000000 & x_0=-0.577350269 & \sim f^{(4)}(x) \\ & w_1=1.000000000 & x_1=+0.577350269 & \\ 3 & w_0=0.555555556 & x_0=-0.774596669 & \sim f^{(6)}(x) \\ & w_1=0.888888889 & x_1=+0.000000000 & \\ & w_2=0.555555556 & x_2=+0.774596669 & \\ 4 & w_0=.0.347854845 & x_0=-0.861136312 & \sim f^{(8)}(x) \\ & w_1=0.652145155 & x_1=-0.339981044 & \\ & w_2=0.652145155 & x_2=+0.339981044 & \\ & w_3=0.347854845 & x_3=+0.861136312 & \\ 5 & w_0=0.236926885 & x_0=-0.906179846 & \sim f^{(10)}(x) \\ & w_1=0.478628670 & x_1=-0.538469310 & \\ & w_2=0.568888889 & x_2=+0.000000000 & \\ & w_3=0.478628670 & x_3=+0.538469310 & \\ & w_4=0.236926885 & x_4=+0.906179846 & \\ 6 & w_0=0.171324492 & x_0=-.0.932469514 & \sim f^{(12)}(x) \\ & w_1=0.360761573 & x_1=-0.661209386 & \\ & w_2=0.467913935 & x_2=-0.238619186 & \\ & w_3=0.467913935 & x_3=+0.238619186 & \\ & w_4=0.360761573 & x_4=+0.661209386 & \\ & w_5=0.171324492 & x_5=+0.932469514 & \\ \hline \end{array} \end{split}\]
\[ I=\int_{-1}^1 f(x) d x \approx I_G=\sum_{i=0}^{n-1} w_i f\left(x_i\right) \]

Пример применения квадратуры Гаусса#

Вычислим интеграл

\[ I=\int_0^\pi 4 x^3 d x=\left.x^4\right|_0 ^\pi=\pi^4=97,409091 \]

Используя 2-х точечную квадратуру Гаусса-Лежандра, мы получим численный ответ $\( I \approx I_G=\frac{(\pi-0)}{2}\left[4\left(\frac{\pi+0}{2}-\frac{\pi-0}{2 \sqrt{3}}\right)^3+4\left(\frac{\pi+0}{2}+\frac{\pi-0}{2 \sqrt{3}}\right)^3\right]=97,409091 \)$

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

Пример на scipy#

Узлы, как и коэффициенты, уже за нас вычислены. Достаточно обратиться к нужным модулям.

import numpy as np
from scipy.special import roots_legendre

roots, weights = roots_legendre(9) # получаем x_i и w_i для нужного n

f = lambda t : t + 1/t # интеграл этой функции на [1, 2] равен 1.5 + np.log(2)

# Преобразуем корни многочленов Лежандра от стандартного [-1, 1] до [a, b]

a = 1
b = 2

t = (b - a)/2 * roots + (a + b)/2

print(np.sum((b - a)/2 * f(t) * weights))
print(1.5 + np.log(2))

2.1931471805599276
2.1931471805599454

Получили ответ с точностью до 14-го знака после запятой, зная значение функции только в 9-ти точках…

Достоинства и недостатки квадратур Гаусса#

Достоинства

  • Требуется малое количество вычислений подынтегральной функции (обычно \(n=10\) хватает для всего)

  • Маленькая погрешность округления из-за малого количества арифметических действий

  • Скорость-скорость-скорость

Недостатки

  • Требуется вычисление функции в определенных точках, что не всегда осуществимо на практике

  • Не подходит для осциллирующих функций

Квадратуры Гаусса-Кристофеля#

В некоторых приложениях возникает необходимость вычисления интегралов с заданной весовой функцией \(w(x)\):

\[\int_a^b f(x) w(x) d x\]

При некотором наборе весовых функций могут оказаться полезными квадратурные формулы Гаусса-Кристоффеля вида

\[ \int_a^b f(x) w(x) d x=\sum_{i=1}^n c_i f\left(x_i\right) . \]

Теорема о точности квадратур Гаусса.

Квадратурная формула выше точна для многочленов степени не выше \(2 n-1\), если

  • ее узлами \(x_i\) являются корни многочлена \(p_n(x)\) из семейства многочленов, ортогональных на промежутке интегрирования \((a, b)\) с весом \(w(x)\):

\[ \left\langle p_k \cdot p_m\right\rangle=\int_a^b p_k(x) p_m(x) \omega(x) d x=0, \quad m \neq k \]

  • весовыми коэффициентами являются числа

\[ c_i=\int_{a}^b \frac{\left(x-x_1\right)\left(x-x_2\right) \ldots\left(x-x_{i-1}\right)\left(x-x_{i+1}\right) \ldots\left(x-x_n\right)}{\left(x_i-x_1\right)\left(x_i-x_2\right) \ldots\left(x_i-x_{i-1}\right)\left(x_i-x_{i+1}\right) \ldots\left(x_i-x_n\right)} w(x) d x . \]

Таблица с видами весового коэффициента, и какие к ним брать многочлены (и на каком отрезке).

\(p_n(x)\)

\(\omega(x)\)

\(a, b\)

Legendre

\(P_n(x)\)

1

-1,1

Hermite

\(H_n(x)\)

\(e^{-x^2}\)

- ထ, ထ

Chebyshev I kind

\(T_n(x)\)

\(\frac{1}{\sqrt{1-x^2}}\)

-1,1

Chebyshev II kind

\(U_n(x)\)

\(\sqrt{1-x^2}\)

-1,1

Laguerre

\(L_n^{(\alpha)}(x)\)

\(x^\alpha e^{-x}\)

0, ထ
Содержание