Квадратурная формула Гаусса: узлы и веса

Квадратурная формула Гаусса - это способ приближённо вычислить определённый интеграл $\int_a^b f(x)\,dx$ суммой $\sum w_i f(x_i)$, в которой и узлы $x_i$, и веса $w_i$ выбраны оптимально, а не заданы заранее. В отличие от формул Ньютона-Котеса (трапеции, Симпсон), где узлы равноотстоящие, в гауссовой квадратуре узлы - это корни ортогональных многочленов, и за счёт свободы их выбора $n$-узловая формула точна для всех многочленов степени до $2n-1$. Ниже разберём, откуда берутся узлы, как считаются веса, что такое алгебраическая точность, как пересчитать формулу с канонического отрезка $[-1,1]$ на произвольный $[a,b]$ и как оценить остаток.

Идея: $2n$ свободных параметров вместо $n$

В формуле Ньютона-Котеса узлы фиксированы (равномерная сетка), свободны только $n$ весов - поэтому такая формула точна для многочленов степени $n-1$ (или $n$ при чётном числе узлов). Гаусс заметил: если разрешить выбирать ещё и узлы, то свободных параметров становится $2n$ ($n$ узлов и $n$ весов). Накладывая $2n$ условий точного интегрирования мономов $1, x, x^2, \ldots, x^{2n-1}$, получаем формулу с алгебраической точностью $2n-1$:

$\int_{-1}^{1} f(x)\,dx \approx \sum_{i=1}^{n} w_i\, f(x_i), \qquad \text{точно при } f \in \mathbb{P}_{2n-1}.$

Это почти вдвое больше, чем у Ньютона-Котеса при том же числе вычислений $f$. Цена - узлы перестают быть «круглыми» числами: они иррациональны и табулируются заранее.

Tool: узлы и веса квадратуры Гаусса

Чтобы не выписывать корни полиномов Лежандра и не решать систему на веса вручную, опишите задачу в форме ниже: число узлов, отрезок интегрирования и (если нужно) подынтегральную функцию. Получите таблицу узлов $x_i$ и весов $w_i$, пересчёт на ваш $[a,b]$, приближённое значение интеграла и оценку погрешности.

Узлы как корни полиномов Лежандра

Ключевой факт: узлы $x_i$ классической квадратуры Гаусса-Лежандра на $[-1,1]$ - это корни многочлена Лежандра $P_n(x)$. Многочлены Лежандра ортогональны на $[-1,1]$ с весом $1$:

$\int_{-1}^{1} P_m(x)\,P_n(x)\,dx = 0, \qquad m \ne n.$

Почему именно корни $P_n$? Любой многочлен $f$ степени $\le 2n-1$ делится на $P_n$ с остатком: $f = q\,P_n + r$, где $\deg q \le n-1$, $\deg r \le n-1$. Интеграл от $q P_n$ равен нулю по ортогональности ($q$ раскладывается по $P_0,\ldots,P_{n-1}$), а в узлах $P_n(x_i)=0$, значит $f(x_i) = r(x_i)$. Остаётся проинтегрировать $r$ - многочлен степени $\le n-1$, который $n$-узловая интерполяционная формула берёт точно. Так свобода выбора узлов «гасит» старшую половину степеней.

Все $n$ корней $P_n$ вещественны, различны и лежат строго внутри $(-1,1)$, симметрично относительно нуля. Для малых $n$ они известны в радикалах:

$n=2:\ x = \pm \tfrac{1}{\sqrt{3}}; \qquad n=3:\ x = 0,\ \pm\sqrt{\tfrac{3}{5}}.$

Веса квадратуры и их свойства

После того как узлы найдены, веса $w_i$ определяются из условия точного интегрирования базисных многочленов. Замкнутая формула для весов Гаусса-Лежандра:

$w_i = \frac{2}{\left(1 - x_i^2\right)\bigl[P_n'(x_i)\bigr]^2}.$

Эквивалентно, $w_i = \int_{-1}^{1} \ell_i(x)\,dx$ - интеграл от базисного многочлена Лагранжа $\ell_i$, построенного по узлам $x_i$. Веса обладают важными свойствами:

• Все веса положительны, $w_i > 0$ - это гарантирует устойчивость формулы (нет роста ошибок округления, как в Ньютоне-Котесе высоких порядков, где появляются отрицательные веса).
• Сумма весов равна длине отрезка: $\sum_{i=1}^n w_i = 2$ для $[-1,1]$ (точное интегрирование $f\equiv 1$).
• Веса симметричны: $w_i = w_{n+1-i}$, как и узлы.

Пересчёт на произвольный отрезок $[a,b]$

Таблицы дают узлы и веса для канонического отрезка $[-1,1]$. Чтобы интегрировать по $[a,b]$, делают линейную замену переменной $x = \frac{b-a}{2}\,t + \frac{a+b}{2}$, переводящую $t\in[-1,1]$ в $x\in[a,b]$. Тогда $dx = \frac{b-a}{2}\,dt$, и формула принимает вид:

$\int_a^b f(x)\,dx \approx \frac{b-a}{2}\sum_{i=1}^{n} w_i\, f\!\left(\frac{b-a}{2}\,t_i + \frac{a+b}{2}\right),$

где $t_i$ и $w_i$ - табличные узлы и веса на $[-1,1]$. Множитель $\frac{b-a}{2}$ - якобиан замены; узлы сдвигаются в центр отрезка и масштабируются. Эта же замена объясняет, почему аппроксимацию полиномами Чебышёва и другие интерполяционные методы тоже строят сначала на $[-1,1]$, а затем переносят на рабочий отрезок.

Алгебраическая точность и оценка погрешности

Главное достоинство формулы - алгебраическая точность $2n-1$: она интегрирует точно любой многочлен этой степени. Для произвольной гладкой $f$ остаток выражается через производную порядка $2n$:

$R_n = \frac{(b-a)^{2n+1}\,(n!)^4}{(2n+1)\bigl[(2n)!\bigr]^3}\, f^{(2n)}(\xi), \qquad \xi \in (a,b).$

Из формулы видно: если $f^{(2n)}$ ограничена, погрешность убывает чрезвычайно быстро с ростом $n$ - у гауссовой квадратуры спектральная сходимость на аналитических функциях, в отличие от полиномиального $O(h^2)$ у трапеций. Сравните с оценкой остатка метода трапеций, где порядок фиксирован и точность набирают только дроблением шага. Для $n=2$ остаток пропорционален $f^{(4)}$, как у Симпсона, но с меньшей константой и при двух узлах вместо трёх.

Когда применять и составная формула

Гауссова квадратура выигрывает, когда вычисление $f$ дорого (каждая точка на счету) и функция гладкая. Если же $f$ имеет особенности, разрывы производных или задана таблично на фиксированной сетке - узлы Гаусса (иррациональные, не совпадающие с узлами таблицы) применить нельзя, и тогда удобнее Ньютон-Котес. На длинном отрезке или при умеренной гладкости используют составную гауссову квадратуру: $[a,b]$ режут на подотрезки и на каждом применяют формулу с малым $n$ (обычно $n=2$ или $3$). Существуют и родственные правила: Гаусса-Лобатто (узлы включают концы $\pm 1$), Гаусса-Радо (один конец), Гаусса-Кронрода (добавочные узлы для апостериорной оценки ошибки) и квадратуры с весовыми функциями - Гаусса-Эрмита, Гаусса-Лагерра, Гаусса-Чебышёва для интегралов с весами $e^{-x^2}$, $e^{-x}$, $(1-x^2)^{-1/2}$.

Частые ошибки

• Забывают якобиан $\frac{b-a}{2}$ при переходе на $[a,b]$. Узлы пересчитывают, а множитель перед суммой опускают - ответ отличается в $\frac{b-a}{2}$ раз. На $[-1,1]$ множитель равен $1$ и незаметен, поэтому ошибка вылезает только на других отрезках.
• Считают, что $n$-узловая формула точна для степени $n-1$. Это про Ньютон-Котес. У Гаусса точность $2n-1$ - вдвое выше, в этом весь смысл выбора узлов.
• Берут равноотстоящие узлы и называют это «формулой Гаусса». Узлы Гаусса - корни $P_n(x)$, они неравномерны и сгущаются к концам. Равномерная сетка даёт обычный Ньютон-Котес.
• Применяют к негладкой функции и удивляются плохой точности. Остаток зависит от $f^{(2n)}$; если высоких производных нет (излом, особенность), спектральная сходимость теряется - нужна составная формула или дробление в точке особенности.
• Путают веса с длинами подотрезков. $w_i$ - не «ширина», а коэффициенты квадратуры; их сумма равна $2$ на $[-1,1]$ (или $b-a$ после пересчёта), но по отдельности они не равны расстояниям между узлами.

FAQ

Почему узлы - именно корни полиномов Лежандра, а не другие точки? Потому что только при таком выборе «лишняя» старшая половина степеней многочлена обнуляется: разложив $f = qP_n + r$, мы используем ортогональность $P_n$ ко всем многочленам меньшей степени (интеграл от $qP_n$ равен нулю) и обращение $P_n$ в ноль в узлах ($f(x_i)=r(x_i)$). Любой сдвиг узлов разрушает это и роняет точность с $2n-1$ до $n-1$.

Как получить узлы и веса для большого $n$, если в радикалах их не выписать? Стандарт - алгоритм Голуба-Уэлша: узлы суть собственные значения симметричной трёхдиагональной матрицы Якоби (её элементы - коэффициенты трёхчленного рекуррентного соотношения для $P_n$), а веса выражаются через первые компоненты собственных векторов. На практике берут готовые таблицы или вызывают библиотечную функцию (например, numpy.polynomial.legendre.leggauss).

Чем формула Гаусса лучше Симпсона при том же числе точек? При $n$ узлах Гаусс точен до степени $2n-1$, а Симпсон (как Ньютон-Котес) - до степени $\sim n$. При двух узлах Гаусс уже даёт точность 3-й степени, тогда как Симпсону нужно три узла для той же 3-й степени. На гладких функциях Гаусс сходится спектрально, а Симпсон - лишь как $O(h^4)$.

Коротко

Квадратурная формула Гаусса приближает $\int_a^b f(x)\,dx$ суммой $\frac{b-a}{2}\sum w_i f(x_i)$, в которой узлы $x_i$ - корни полинома Лежандра $P_n$ на $[-1,1]$, а веса $w_i = \frac{2}{(1-x_i^2)[P_n'(x_i)]^2}$ положительны и в сумме дают $2$. Свобода выбора и узлов, и весов даёт алгебраическую точность $2n-1$ - вдвое выше Ньютона-Котеса при том же числе вычислений функции; остаток пропорционален $f^{(2n)}(\xi)$, что обеспечивает спектральную сходимость на гладких функциях. Для отрезка $[a,b]$ применяют линейную замену $x = \frac{b-a}{2}t + \frac{a+b}{2}$ с якобианом $\frac{b-a}{2}$. На негладких или таблично заданных функциях метод теряет преимущество - там используют составную формулу или правила Ньютона-Котеса.