Численное интегрирование

Значение определенного интеграла $$ \int_a^b f(x)\, dx $$ можно приближенно вычислить, используя одну из трех квадратурных формул: формулу прямоугольников, формулу трапеций, формулу Симпсона. Во всех них мы делим отрезок $[a, b]$ на $n$ равных частей и вычисляем интеграл, используя только значения функции в точках $$ x_n = a + i\cdot h, \quad i=0, 1, 2, \dots, n $$ Через $h$ мы обозначили длину отрезков разбиения: $$ h = (b-a)/n $$

Наиболее проста формула прямоугольников: $$ \int_a^b f(x)\, dx \approx h\sum_{i=0}^{n-1}f(x_i) $$ Однако точность вычисления интеграла по этой формуле невысока. Формулу прямоугольников имеет смысл использовать лишь для периодических функций при интегрировании по периоду.

Заметно более высокую точность обеспечивает формула трапеций, при том, что скорость вычисления по ней такая же, да и ее сложность ничуть не выше, чем у формулы прямоугольников: $$ \int_a^b f(x)\, dx \approx h(\frac{1}{2}(f(a)+f(b)) + \sum_{i=1}^{n-1}f(x_i)) $$ Основана она на том, что интеграл по каждому маленькому отрезку разбиения можно заменить на площадь трапеции: $$ \int_{x_i}^{x_{i+1}} f(t)\, dt \approx \frac{1}{2}h(f(x_i)+f(x_{i+1})) $$ Формула трапеций получается суммированием по всем отрезкам разбиения.

Отметим, что формула прямоугольников дает точное значение интеграла для функций-констант, а формула трапеций — для линейных функций. Давайте выведем квадратурную формулу, которая вычисляет точное значение интеграла для полиномов второй степени (квадратичных функций). Эта формула называется Формулой Симпсона. Как окажется, она точно не только для квадратичнеых функций, но даже и для полиномов третьей степени.

Формула Симпсона основана на следующей простой формуле. Рассмотрим отрезок $[x_0, x_2]$ и обозначим через $x_1$ середину этого отрезка: $x_1 = (x_0 + x_2)/2$. Обозначим $h = x_1 - x_0$. Пусть функция $f(x)$ принимает в точках $x_0$, $x_1$, $x_2$ значения $y_0$, $y_1$, $y_2$: $$ f(x_i) = y_i, \quad i = 0, 1, 2. $$

Приблизим функцию $y = f(x)$ на отрезке $[x_0, x_2]$ параболой $y = g(x) = \alpha x^2 + \beta x + \gamma$, принимающей в точках $x_i$ значения $y_i$. Тогда площадь под графиком параболы равна $$ \int_{x_0}^{x_2} g(x) dx = \frac{1}{6}(x_2 - x_0) (y_0 + 4y_1 + y_2) = \\ \frac{1}{3}h(y_0 + 4y_1 + y_2) $$

Эта формула легко запоминается: значения функции на концах отрезка суммируются с весовыми коэффициентами 1, значение в середине отрезка — с весовым коэффициентом 4, сумма делится на 6 (это сумма весовых коэффициентов) и умножается на длину отрезка интегрирования. Назовем ее условно малой формулой Симпсона. Ее нетрудно вывести, выписав коэффициенты квадратного трехчлена с помощью интерполяционной формулы Лагранжа и вычислив интеграл от полученного многочлена. Выкладки, однако, получаются довольно громоздкими, поэтому удобнее доказать данную формулу следующим образом. Заметим, что результат вычисления по этой формуле линено зависит от функции $g(x)$, поэтому достаточно доказать ее для степеней переменной: $g(x) = x^k$, $k = 0, 1, 2$. Для $k = 0, 1$ формула очевидна. Проверим для $k = 2$. Вычислим сначала интеграл: $$ \int_{x_0}^{x_2} x^2 dx = \int_{-h}^{h} (x-x_1)^2\, dx = \frac{1}{3} (x-x_1)^3 \bigg|_{-h}^{h} = \\ \frac{1}{3}\big( (h - x_1)^3 - (-h - x_1)^3 \big) = \frac{1}{3}\big( (h - x_1)^3 + (h + x_1)^3 \big) = \\ \frac{1}{3} (2h^3 + 6hx_1^2) $$ А теперь применим малую формулу Симпсона для функции $y=x^2$: $$ \frac{1}{3}h (y_0 + 4y_1 + y_2) = \frac{1}{3}h \big((x_1-h)^2 + 4x_1^2 + (x_1+h)^2 \big) = \\ \frac{1}{3}h (6x_1^2 + 2h^2) = \frac{1}{3} (6hx_1^2 + 2h^3) $$ Получили тот же результат. Итак, малая формула Симсона доказана для многочленов степени не выше 2. Удивительно, однако, что эта формула дает точное значение определенного интеграла и для многочленов степени 3! (Это следствие симметричности используемой весовой суммы.) Доказать это утверждение можно аналогично, проверив ее для функции $y = x^3$ (здесь мы эту проверку опускаем).