можно ли вычислить интеграл с заданной погрешностью

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

Оглавление

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

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

Пусть требуется вычислить определенный интеграл

Суть методов численного интегрирования состоит в замене подынтегральной функции f(x) вспомогательной, интеграл от которой легко вычисляется в элементарных функциях. Обычно f(x) заменяется некоторым интерполяционным многочленом, что приводит к квадратурным формулам:

,

где x_i – узлы интерполяции;

i – произвольный номер узла;

R – остаточный член, или погрешность метода.

Неучет (отбрасывание) R приводит к погрешности усечения. К этим погрешностям в процессе вычислений добавляются погрешности округления. Геометрическая интерпретация численного интегрирования представлена на рис. 5.1. и 5.2.

1. Методы прямоугольников

Самыми простыми методами численного интегрирования являются методы прямоугольников. При этом непосредственно используется замена определенного интеграла интегральной суммой:

, h_i = x_i – x_{i —1} , x_{i —1} ≤ z_i ≤ x_i

;

Для частного случая h_i = h = const формулы примут вид

;

» width=»147″ height=»60″ alt=» » />;

Если координату выразить через начальную точку и принять, что I ≈ S, то получим формулы, готовые для применения в операторе цикла с переменной:

» width=»188″ height=»60″ alt=» » /> – для метода левых прямоугольников;

» width=»153″ height=»60″ alt=» » /> – для метода правых прямоугольников;

» width=»203″ height=»52″ alt=» » /> – для метода средних прямоугольников.

» width=»480″ height=»308.08510638298″ alt=» » style=»display: block; margin-left: auto; margin-right: auto;» />

Рис. 5.1. Графическая интерпретация методов прямоугольников

2. Метод трапеций

После сложения этих соотношений получим формулу трапеций

» width=»212″ height=»52″ alt=» » />.

Если шаг интегрирования постоянный (h_i = h = const), то

» width=»241″ height=»60″ alt=» » />.

Если координату выразить через начальную точку и принять, что I » S, то получим формулу, готовую для применения в операторе цикла с переменной

» width=»325″ height=»52″ alt=» » />.

» width=»480″ height=»308.33333333333″ alt=» » />

Рис. 5.2. Графическая интерпретация метода трапеций

3. Процедура вычисления интеграла

Процедура вычисления интеграла при заданном числе разбиений и границами интегрирования должна иметь входные параметры этих величин и выходной параметр результата. Выходной параметр должен содержать перед собой ключевое слово var для передачи параметра по ссылке.

Здесь должен вычисляться шаг в соответствии с числом разбиений.

Сумма может вычисляться при помощи цикла for по числу разбиений или при помощи циклов repeat или while с перемещением текущей точки от начальной до конечной. Например, для метода правых прямоугольников процедура может быть описана так:

procedure integral(a,b:real;n:integer; var s:real);

var i:integer; h:real;

begin

h:=abs(b-a)/n; s:=0;

for i:=1 to n do

s:=s+h*f(a+h*i);

Конечно, перед описанием процедуры необходимо описать подынтегральную функцию f(x). Смотрите описание функций в 3.1.

Замечание! При вычислении интеграла функция должна браться по модулю, т. к. площадь фигуры всегда положительна, независимо от того, что функция может быть и отрицательной. Для упрощения программы модуль можно вставить сразу в описание функции.

4. Вычисление интегралов с заданной точностью и оценка методов интегрирования

Оценка погрешности усечения R в формулах численного интегрирования оказывается трудоемкой и малоэффективной из-за трудностей оценки производных высокого порядка подынтегральных функций. Поэтому в практических расчетах для достижения требуемой точности вычислений или допустимой погрешности E используют правило Рунге.

Согласно этому правилу вычисление заданного интеграла проводят для разных интервалов разбиения отрезка [A,B]. Так, если начальное число интервалов разбиения есть n и соответствующее ему приближенное значение интеграла – I_n, то для числа интервалов 2n получим значение интеграла I_2n. Число интервалов можно увеличивать в 2, 3 и т. д. раз по сравнению с базовым значением n. При двукратном увеличении числа отрезков погрешность D приближенного значения интеграла для методов прямоугольников и трапеций оценивается как

Если Δ > E, то количество интервалов разбиения [A,B] опять увеличивают вдвое, т. е. значение I вычисляют для 4n. Такое удвоение повторяют до тех пор, пока не выполнится условие Δ Основная часть программы должна содержать ввод с клавиатуры начального числа разбиений и погрешности вычислений. Для обеспечения заданной точности интеграла процедура его вычисления должна выполняться в цикле. Причём сначала сохраняется старое значение интеграла, затем увеличивается число разбиений и вычисляется новое значение интеграла. Цикл работает пока разница между старым и новым значениями интеграла не станет меньше заданной точности.

При выводе результата полезно также вывести конечное число разбиений, чтобы знать, на сколько частей нужно разбить отрезок для обеспечения заданной точности вычисления интеграла.

Источник

Оценка погрешности численного интегрирования

Тема 1.4. Численное интегрирование

1.4.1. Постановка задачи

1.4.2. Метод прямоугольников

1.4.3. Формула трапеций

1.4.4. Формула Симпсона

1.4.5. Оценка погрешности численного интегрирования

1.4.6. Тестовые задания по теме «Численное интегрирование»

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

Из курса математического анализа известно, что, если функция f(x) непрерывна на отрезке [a;b] и дифференцируема, то определенный интеграл от этой функции в пределах от a до b существует и может быть вычислен по формуле Ньютона-Лейбница:

Если первообразную функцию F(x) не удается выразить аналитически через элементарные функции или если при проведении практических расчетов подынтегральная функция f(x) задается в виде таблицы, то это приводит к необходимости замены аналитического интегрирования численными методами.

Для функции f(x), заданной в прямоугольной системе координат на интервале [a;b], этот интеграл численно равен площади, ограниченной кривой f(x), осью Ox и двумя ординатами ac и bd.

Применительно к однократному интегралу, формулы численного интегрирования представляют собой квадратурные формулы вида:

Искомый определенный интеграл можно представить в виде суммы интегралов:

На каждом i-м отрезке функция аппроксимируется (заменяется) некоторой другой легко интегрируемой функцией g_i(x). В результате получаем следующую квадратурную формулу:

.

Для решения поставленной задачи подынтегральную функцию f(x) необходимо заменить приближенной функцией, которая может быть проинтегрирована в аналитическим виде. В качестве такой функции обычно используют полином Р(х) с узлами интерполяции в точках х₀, х₁, х₂, …,х_n. В этих точках значения функции и интерполяционного полинома полностью совпадают f(x_i) = Р(x_i).

Для получения простых формул интегрирования используют полиномы нулевой, первой и второй степени и соответственно получают формулы численного интегрирования: прямоугольников, трапецийиСимпсона.

Очевидно, что замена функции f(x) интерполирующим полиномом приводит к образованию погрешности вычисления значения интеграла

где I₁ – точное значение интеграла, I – значение интеграла, вычисленного численным методом, а – погрешность метода.

Отметим, что увеличение числа подынтервалов n (или уменьшение длины шага интегрирования h) ведет к уменьшению погрешности.

Метод прямоугольников

Заменим подынтегральную функцию f(x) в пределах элементарного отрезка [x_i;x_i+1] интерполяционным многочленом нулевой степени (рис.1.4.2-1), то есть постоянной величиной, равной либо f(x_i), либо f(x_i+1).

(1.4.2-1)
(1.4.2-2)

Формула (1.4.2-1) называется формулой левых прямоугольников, а формула
(1.4.-2.2) – формулой правых прямоугольников.

Для вычисления определенного интеграла может быть использована и формула средних прямоугольников (1.4.2-3), в которой на элементарном отрезке интегрирования функция f(x)тоже заменяется интерполяционным многочленом нулевой степени, но равным значению функции в середине отрезка:

(1.4.2-3)

Схема алгоритма метода средних прямоугольников приведена на рис. 1.4.2-2.

Рис. 1.4.2-2. Схема алгоритма интегрирования по методу средних прямоугольников с

использованием правила Рунге

Формула трапеций

Разобьем интервал интегрирования [a;b] на n равных отрезков (рис. 1.4.3-1) и восстановим из полученных точек a, х₁, x₂, …, b перпендикуляры до пересечения с графиком функции. Соединив последовательно точки пересечения, представим площадь полученной криволинейной трапеции как сумму прямолинейных трапеций, площади которых легко подсчитать. Заменив подынтегральную функцию f(x) в пределах элементарного отрезка [x_i;x_i+1] интерполяционным многочленом первой степени, получим следующие формулы для элементарных площадей:

Тогда общая площадь равна:

Отсюда получаем формулу трапеций:

(1.4.3-1)

Схема алгоритма метода трапеций приведена на рис. 1.4.3-2.

Рис. 1.4.3-2. Схема алгоритма интегрирования по методу трапеции с использованием

Формула Симпсона

Для получения формулы Симпсона применяется квадратичный интерполирующий полином, следовательно, за элементарный интервал интегрирования принимается отрезок [x_i;x_i+2]. Поэтому разобьем интервал интегрирования [a;b] наn отрезков, где n=2m – четное число (рис. 1.4.4-1).

Для получения интерполирующей функции на интервале [x_i;x_i+2] воспользуемся первой интерполяционной формулой Ньютона, используя в качестве узлов интерполяции точки x_i, х_i+1 и x_i+2.

(1.4.4-1)

В пределах отрезка [x_i;x_i+2], на котором подынтегральная функция аппроксимирована многочленом (1.4.4-1), получим приближенную формулу Симпсона:

(1.4.4-2)

Для отрезка [x₀;x₂]

Для отрезка [x₂;x₄]

Тогда для всего интервала интегрирования [a;b] формула Симпсона выглядит следующим образом:

(1.4.4-3)

при

Схема алгоритма метода Симпсона приведена на рис. 1.4.4-2.

Рис. 1.4.4-2. Схема алгоритма интегрирования по методу Симпсона с использованием

Оценка погрешности численного интегрирования

Замена подынтегральной функции интерполяционным полиномом приводит к погрешности вычисления его значения R=|I₁ – I|, где

Очевидно, что вычислить эту погрешность можно только, если известно точное значение интеграла. Поэтому на практике принято проводить оценку погрешности численного интегрирования следующим образом (подынтегральная функция задана таблично (Т) или аналитически (А)):

· при использованииформул левых или правых прямоугольников

· при использованииформулы средних прямоугольников

· при использованииформулы трапеций

· при численном интегрировании по формуле Симпсона:

В приведенных выше формулах: a, b–границы интервала интегрирования; h=(b-a)/n–шаг интегрирования; . и – среднее арифметическое, соответственно, первых, вторых и четвертых конечных разностей.

Поскольку в формуле погрешности для метода трапеций присутствует вторая производная, а в формуле Симпсона – четвертая, то формула трапеций точна только для линейных функций, а формула Симпсона для линейных, квадратичных и кубических.

Из приведенных формул видно, что уменьшение шага интегрирования (h) приводит к уменьшению погрешности. При этом, поскольку квадратичная интерполяция представляет функцию с большей точностью, чем линейная, то при использовании формулы Симпсона требуемая точность достигается при меньших значениях n(количестве разбиений), чем, например, при использовании формулы трапеций и формулы прямоугольников.

Формулы для оценки погрешности могут быть также использованы для выбора числа разбиений n или шага интегрирования h, необходимых для обеспечения заданной точности. Однако практическое использование этих формул ограничено в связи с трудоемкостью их вычислений, поэтому при реализации численных методов на ПК используется прием, позволяющий получить оценку погрешности в неявном виде. Этот прием основан на двукратном вычислении значения интеграла вначале с шагом h(где h=(b-a)/n), а затем с шагом h/2. Полученные значения интегралов I_h и I_h/2 могут быть применены для оценки погрешности интегрирования по формуле:

(1.4.5-1)

где: k=1 – для формул левых и правых прямоугольников;

k=2 – для формул трапеции и средних прямоугольников;

k=4 – для формулы Симпсона.

Если полученная погрешность не удовлетворяет требуемой точности, то вычисляется значение интеграла при h=h/4 и снова оценивается погрешность, и т.д. до тех пор, пока не окажется, что погрешность стала меньше заданной точности. Это правило называется правилом Рунге (или правилом двойного просчета).

Пример 1.4.5-1. Вычислить значение определенного интеграла

Предположим, что, подынтегральная функция задана таблично:

Используем формулы правыхилевых прямоугольников, считая, что h = 0.2, а n=(b-a)/h=5, имеем:

Используем формулу трапеций:

Используем формулы средних прямоугольников:

Используем формулу Симпсона приm=2∙n=10(2∙5) и шаге h=0.1:

Произведем оценку погрешности каждого из полученных значений, используя известное аналитическое выражение подынтегральной функции f(x):

Следовательно,

Анализируя значения погрешностей, можно с уверенностью сказать, что самый точный результат получен с использованием формулы Симпсона.

1.4.6. Тестовые задания по теме
«Численное интегрирование»

1. Численное значение интеграла равно

1)площади, ограниченной кривой f(x), осью 0x и двумя ординатами в точках a и b

2)площади прямоугольника

3)площади прямоугольной трапеции

4)в списке нет правильного ответа

Источник

C++ и Численные Методы: Приближенное интегрование по Ньютону-Котесу

Методы Ньютона-Котеса — это совокупность техник приближенного интегрирования, основанных на:

Метод трапеций

Метод трапеций — простейший из рассмотренных. В качестве примера возьмем следующий интеграл:

Точность приближения зависит от числа N отрезков, на которые разбивается промежуток интегрирования. Таким образом, длина промежутка:

Площадь трапеции может быть вычислена по формуле:

Суммируя все вышесказанное, приближенное значение интеграла вычисляется по формуле:

Функция, вычисляющая интеграл методом трапеций должна принимать 4 параметра:

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

Метод Симпсона заключается в интегрировании интерполяционного многочлена второй степени функции f(x) с узлами интерполяции a, b и m = (a+b)/2 — параболы p(x).Для повышения точности имеет смысл разбить отрезок интегрирования на N равных промежутков(по аналогии с методом трапеций), на каждом из которых применить метод Симпсона.

Площадь параболы может быть найдена суммированием площадей 6 прямоугольников равной ширины. Высота первого из них должна быть равна f(a), с третьего по пятый — f(m), шестого — f(m). Таким образом, приближение методом Симпсона находим по формуле:

Метод Ромберга

Пусть T(x) — приближение интеграла, полученное методом трапеций с шагом x. Получим 3 таких приближения, уменьшая размер шага в 2 раза при каждом вычислении.

Построим теперь симметричную относительно оси y параболу, проходящую через точки T(1) и T(1/2) чтоб экстраполировать полученные значения для x стремящегося к 0.

Следовательно, каждый член первого столбца R(n, 0) приближений Ромберга эквивалентен решениям полученным методом трапеций, а каждое решение второго столбца R (n, 1) — методом Симпсона. Таким образом, формулы для приближенного интегрирования методом Ромберга:

Источник