Метод Симпсона: численное интегрирование

Метод Симпсона - квадратурная формула численного интегрирования, в которой подынтегральная функция на каждом маленьком отрезке заменяется параболой через три равноотстоящие точки, а интеграл от параболы берётся аналитически. В обмен на чуть более громоздкую формулу, чем у трапеций, метод получает алгебраическую точность 3 (точно интегрирует все многочлены до третьей степени) и глобальную погрешность $O(h^4)$ - на порядок лучше трапеций при том же числе вычислений $f$.

Постановка задачи и идея

Нужно вычислить определённый интеграл $I = \int_a^b f(x)\, dx$, когда первообразная неизвестна, выражается через специальные функции (как $\int e^{-x^2} dx$) или функция вообще задана таблично. Идея всех квадратурных формул одна: разбить $[a, b]$ на куски, на каждом заменить $f$ простой функцией и сложить. Прямоугольники заменяют $f$ константой, трапеции - линейной, метод Симпсона - параболой $P_2(x)$ через три точки.

Возьмём отрезок $[a, b]$ и три узла: концы $x_0 = a$, $x_2 = b$ и середину $x_1 = (a+b)/2$. Шаг $h = (b - a)/2$. Через $f_0, f_1, f_2$ проходит единственный многочлен Лагранжа $P_2(x)$. Заменяем $f$ на $P_2$ и интегрируем:

$\int_a^b P_2(x)\, dx = \frac{h}{3}\bigl[f_0 + 4 f_1 + f_2\bigr].$

Это базовая формула Симпсона для «двойного» отрезка длиной $2h$. Веса $1, 4, 1$ получаются прямым вычислением интегралов от базисных многочленов Лагранжа. Алгебраическая точность - 3: формула точна не только для $1, x, x^2$, но и для $x^3$ - за счёт симметрии узлов относительно $x_1$.

Tool: посчитать интеграл методом Симпсона

Подставлять «руками» формулу $\frac{h}{3}[f_0 + 4f_1 + f_2]$ для $n = 16$ узлов утомительно: восемь четвёрок, семь двоек, два краевых значения. Ниже - мини-форма: выбираешь функцию из списка (или вписываешь свою), границы $a, b$ и число отрезков $n$ - получаешь подробный разбор с подстановкой в составную формулу и сравнением с аналитическим ответом, где он есть.

Составная формула Симпсона

Один глобальный шаг $h = (b - a)/2$ - это грубо. Поэтому разбивают $[a, b]$ на $n$ равных отрезков длиной $h = (b - a)/n$, причём $n$ обязательно чётное, и применяют базовую формулу к каждой паре соседних отрезков. Узлы: $x_i = a + i h$, $i = 0, 1, \ldots, n$, значения $f_i = f(x_i)$. После приведения подобных получается составная формула Симпсона:

$I \approx \frac{h}{3}\bigl[f_0 + 4(f_1 + f_3 + \ldots + f_{n-1}) + 2(f_2 + f_4 + \ldots + f_{n-2}) + f_n\bigr].$

Запоминается просто: крайние значения с весом 1, нечётные узлы (центры «двойных» отрезков) - с весом 4, чётные внутренние - с весом 2. В сумме $f_0 + 4\sum_{\text{нечёт}} + 2\sum_{\text{чёт}} + f_n$ участвует $n + 1$ значение функции - столько же, сколько в формуле трапеций, но точность выше в $h^2/h^4 = 1/h^2$ раз.

Оценка погрешности и порядок $O(h^4)$

Точное представление остатка для составной формулы:

$R = -\frac{(b - a)\, h^4}{180}\, f^{(4)}(\xi), \qquad \xi \in (a, b).$

Отсюда оценка сверху:

$|R| \le \frac{(b - a)\, h^4}{180}\, \max_{x \in [a,b]} |f^{(4)}(x)|.$

Главный вывод: погрешность падает как $h^4$. Удвоение числа отрезков (уменьшение $h$ вдвое) снижает ошибку в $2^4 = 16$ раз. У трапеций тот же двукратный шаг даёт 4 раза. Метод Симпсона точен для кубических многочленов: $f^{(4)} \equiv 0$ для любого $f$ степени $\le 3$, остаток обнуляется.

Практически погрешность контролируют правилом Рунге: считают интеграл с шагом $h$ и с шагом $h/2$, разность $I_{h/2} - I_h$ делят на $2^4 - 1 = 15$ - это оценка остатка для $I_{h/2}$. Если разность меньше $\varepsilon$, точность достигнута, а ответ можно уточнить по Ричардсону: $I \approx I_{h/2} + (I_{h/2} - I_h)/15$ - это уже $O(h^6)$. Альтернативный способ связать сумму значений функции с интегралом и явно выписать поправочные члены через производные на концах даёт формула Эйлера-Маклорена - её разложение по сути объясняет, почему симметричные квадратуры типа Симпсона теряют чётные степени $h$ в остатке.

Сравнение с прямоугольниками и трапециями

Все три классические формулы Ньютона-Котеса используют равноотстоящие узлы и отличаются только порядком интерполяции:

Для гладких функций метод Симпсона почти всегда оптимален: точность $h^4$ при тех же $n + 1$ вычислениях $f$, что и у трапеций. Для функций с разрывом производной порядок $h^4$ не достигается - нужна адаптивная схема или замена переменной.

Метод Симпсона 3/8

Если интерполировать тремя интервалами (четырьмя точками $x_0, x_1, x_2, x_3$) кубической параболой и взять интеграл, получится правило Симпсона 3/8:

$\int_{x_0}^{x_3} f\, dx \approx \frac{3h}{8}\bigl[f_0 + 3 f_1 + 3 f_2 + f_3\bigr], \qquad h = (x_3 - x_0)/3.$

Алгебраическая точность та же - 3, порядок тот же - $O(h^4)$, но константа в остатке чуть хуже, и формула требует, чтобы число отрезков было кратно 3. На практике 3/8 идёт добивкой, когда $n$ оказалось нечётным: последние три отрезка считают по 3/8, остальные - обычным Симпсоном.

Адаптивная квадратура Симпсона

Когда $f$ ведёт себя неравномерно - гладко в одной части отрезка и резко меняется в другой, - равномерное разбиение тратит точки впустую. Адаптивный Симпсон работает рекурсивно: для отрезка $[a, b]$ считает $S(a, b)$ на $n = 2$, потом $S(a, m) + S(m, b)$ на середине $m = (a + b)/2$. Если $|S(a, m) + S(m, b) - S(a, b)| < 15 \varepsilon$, считает результат достаточно точным; иначе рекурсивно вызывает себя на половинках с допуском $\varepsilon / 2$. Алгоритм автоматически кладёт больше точек туда, где $|f^{(4)}|$ велика. Это базовый трюк библиотек численного интегрирования (scipy.integrate.quad).

Классические примеры и типовые задачи

Площадь под $\sin x$ на $[0, \pi]$. Точное значение - $\int_0^\pi \sin x\, dx = 2$. Возьмём $n = 4$, $h = \pi/4$: $f_0 = 0$, $f_1 = \sin(\pi/4) \approx 0{,}7071$, $f_2 = 1$, $f_3 \approx 0{,}7071$, $f_4 = 0$. Подстановка: $I \approx \frac{\pi/4}{3}[0 + 4 \cdot 0{,}7071 + 2 \cdot 1 + 4 \cdot 0{,}7071 + 0] \approx 2{,}0046$. Ошибка $0{,}0046$ при всего пяти значениях функции.

Интеграл Гаусса $\int_{-1}^1 e^{-x^2}\, dx$ аналитически выражается через функцию ошибок $\sqrt{\pi}\, \mathrm{erf}(1) \approx 1{,}4936$. Метод Симпсона на $n = 8$ даёт $\approx 1{,}4937$ - никакая первообразная не нужна. Длина дуги $L = \int_a^b \sqrt{1 + (f')^2}\, dx$ почти никогда не берётся аналитически; численный Симпсон на $n = 16$-$32$ даёт точность 4-5 знаков.

В типовых задачах на коллоквиуме просят: посчитать интеграл с заданным $n$ и сравнить с точным значением; оценить минимальное $n$, при котором погрешность $\le \varepsilon$, через формулу $|R| \le \frac{(b-a) h^4}{180} M_4$; применить правило Рунге и уточнить ответ по Ричардсону; сравнить точность Симпсона, трапеций и средних прямоугольников на одном $f$.

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

• Берут нечётное $n$. Составная формула Симпсона не определена для нечётного числа отрезков - нужна либо чётность, либо добивка правилом 3/8 на последних трёх отрезках.
• Путают веса 4 и 2. Узлы с нечётными индексами (середины «двойных» отрезков) - это 4; чётные внутренние - это 2. Проверка: при $f \equiv 1$ сумма должна дать $b - a$.
• Считают шаг через $(b - a)/(n - 1)$. На самом деле $h = (b - a)/n$, а узлов получается $n + 1$. Ошибка сдвигает все $f_i$.
• Применяют Симпсон к функции с особенностью. В точке излома или на $\int_0^1 \sqrt{x}\, dx$ четвёртая производная неограничена, и обещанный $O(h^4)$ не работает - нужно разбить интеграл на гладкие куски или взять адаптивную схему.

FAQ

Почему метод Симпсона точен для кубических многочленов, если интерполируем параболой? Из-за симметрии узлов. Базовая формула берёт три симметричные относительно $x_1$ точки. Кубический член $x^3$ - нечётная функция относительно середины, и его интеграл по симметричному отрезку, посчитанный через симметричные веса, совпадает с точным значением. Формально это видно из остатка $R \propto f^{(4)}(\xi)$: для кубики $f^{(4)} = 0$, и $R = 0$ автоматически.

Когда лучше использовать правило 3/8? Когда число отрезков $n$ кратно трём, но не кратно двум (например, $n = 9, 15, 21$) - тогда обычный Симпсон не применим, а 3/8 ложится ровно. Чаще 3/8 используют как добивку: для $n = 5$ берут $n - 3 = 2$ отрезка обычным Симпсоном плюс последние 3 отрезка - правилом 3/8. На гладких функциях обе формулы дают сопоставимую точность.

Чем метод Симпсона отличается от квадратур Гаусса? Симпсон - это формула Ньютона-Котеса с равноотстоящими узлами. Квадратуры Гаусса выбирают узлы оптимально (нули полиномов Лежандра) и для $k$ узлов дают алгебраическую точность $2k - 1$ против $k - 1$ у Ньютона-Котеса. Гаусс точнее, но узлы и веса заранее не очевидны - нужны таблицы; зато для аналитических функций сходимость экспоненциальная, а не степенная. Симпсон проще объяснить и удобнее, когда $f$ задана таблично с равномерным шагом.

Коротко

Метод Симпсона - квадратурная формула, заменяющая подынтегральную функцию параболой через три равноотстоящие точки. Базовая формула $\frac{h}{3}[f_0 + 4 f_1 + f_2]$ точна для всех многочленов до куба включительно. Составная формула $\frac{h}{3}[f_0 + 4 \sum_{\text{нечёт}} + 2 \sum_{\text{чёт}} + f_n]$ на $n$ чётных отрезках даёт погрешность $O(h^4)$ - на порядок лучше трапеций при том же числе вычислений $f$. Удвоение $n$ снижает ошибку в 16 раз; точность контролируют правилом Рунге, а адаптивные варианты Симпсона лежат в основе библиотек численного интегрирования. Для нечётного числа отрезков работает правило Симпсона 3/8 с тем же порядком. Метод применим к любой гладкой $f$, не требует первообразной и одинаково удобен как для аналитически заданной функции, так и для табличных данных.