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

Формулы численного интегрирования.

При решении многих задач, встречающихся в геометрии, технике, экономике, приходится вычислять определенные интегралы.

Если для подынтегральной функции f (x ) найдена первообразнаяF (x ) , то интеграл, как известно, можно вычислить по формуле Ньютона-Лейбница:

(1)

Однако на практике часто не бывает возможности использовать формулу (1), например, в следующих случаях:

    если первообразная функция F (x ) не выражается в конечном виде через элементарные функции. Это относится, например, к интегралам:

    если аналитическое выражение первообразной функции F (x ) является настолько сложным, что применение формулы (1) становится затруднительным;

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

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

Формула прямоугольников. Вывод этой формулы основан на замене определенного интеграла интегральной суммой. Из анализа известно, что

где
- интегральная сумма для функцииf (x ) на отрезке[ a , b ].

ξ – внутренняя точка отрезка[ a , b ].

Если отрезок [ a , b ] разбить наn равных частей:

а=х 0 , х 1 , …, х п = b ,

х i = = h .

Число h называетсяшагом квадратурной формулы. При этом условии получаем:

Если взять в качестве точек ξ i левые концы частичных отрезков:

f(ξ i ) = f(х i ) (i = 0, 1, …, n-1),

Обозначим f i ) = у i . Заменяя интеграл интегральной суммой, получим приближенное равенство:

, (2)

называемое формулой прямоугольников (с левыми ординатами).

Если взять в качестве точек ξ i правые концы частичных отрезков:

f (ξ i ) = f i ) (i = 1, 2,…, n ),

то получим приближенное равенство:

, (3)

называемое формулой прямоугольников (с правыми ординатами).

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

Пример. Вычислим интеграл, разбив интервал интегрирования на 10 равных частей (n = 10 ). Найдем и запишем в таблицу значения подынтегральной функции

у = в точках деления:

i

х i

у i =

i

х i

у i =

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

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

Значение, полученное по формуле (1):

Мы видим, что формулы прямоугольников дают грубые приближения.

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

Абсолютную погрешность r формул прямоугольников (2) и (3) можно оценить по формуле:

(4)

Идея вывода квадратурных формул трапеций и Симпсона:

подынтегральной функции f ( x ) поставить в соответствие близкую ей функциюg n ( x ) , которую можно проинтегрировать, и приближенно заменить искомый интегралI интегралом от этой функции.

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

Обозначим a = x 0 , b = x 1 .

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

f (x ) у 0 +t у 0 ,

у 0 =f (x 0 ) ,у 1 =f (x 1 ) , у 0 =у 1 - у 0 .

В этом случае

, (5)

Известно, что t =

Отсюда х=х 0 + th и dx =hdt .

При х = х 0 t = 0;

при х =х 1 t = 1 .

Переходя к новой переменной t , получим:

(6)

так как у 0 =у 1 у 0

Формула (6) называется формулой трапеций.

Ее геометрический смысл состоит в том, что на отрезке [х 0 ;х 1 ] криваяу =f(х) заменяется отрезком прямой (хордой), т. е. криволинейная трапеция заменяется прямолинейной.

Значение интеграла, вычисленное по формуле (6), будет равно площади трапеции. На рисунке эта площадь заштрихована.

Как показывает вычислительная практика, при недостаточно малой длине отрезка интегрирования точность результатов, полученных с помощью формулы (6), бывает недостаточной.

Для получения более точного результата поступают следующим образом:

Отрезок интегрирования [а; b ] разбивают на п равных частей точками: х 0 = а, х 1 , х 2 ,…,х n = b . И аппроксимируют кусочно-линейной функцией g п (x) . Применяя формулу (6) на каждом из частичных отрезков интегрирования, получают:

(7)

Сложив равенства, получают формулу, называемую обобщенной формулой трапеций:

(8)

где у i =f(х i ) (i = 0, 1, …, n).

Геометрический смысл этой формулы состоит в том, что кривая - график функции у = f (х) - заменяется ломаной, вписанной в кривую АВ . Площадь криволинейной трапеции заменяется суммой площадей прямолинейных трапеций. Как показывает практика, формула (8) при большом числе точек деления позволяет получать хорошие результаты.

Пример 1. Вычислим по формуле трапеций (8) интеграл , разбив отрезок интегрирования на десять равных частей.

Воспользовавшись данными, занесенными в предыдущую таблицу, получим:

Сравнение полученного результата со значением ln2  0,693147 показывает, что погрешность значения интеграла, вычисленного по обобщенной формуле трапеций, значительно меньше погрешности, допущенной при вычислении этого же интеграла по формуле прямоугольников.

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

(9)

где а < < b ,

а абсолютная погрешность оценивается следующим образом:

(10)

(11)

Формула Симпсона (формула парабол)

Для вычисления интеграла
разобьем отрезок интегрирования на два равных отрезка:

0 , х 1 ] и 1 , х 2 ] (х 0 = а, х 2 =b )

и заменим подынтегральную функцию по формуле квадратичного интерполирования

(12)

где t = .

.

Перейдем к новой переменной интегрирования, учитывая, что

х = х 0 + ht , dx = hdt ,

при х=х 0 t =0

при х=х 2 t =2

(13)

Формула (13) называется формулой Симпсона или формулой парабол.

Ее геометрический смысл состоит в следующем: на отрезке 0 , х 2 ] кривая у = f (x ) заменяется квадратной параболой - графиком интерполяционного многочлена. При вычислении по формуле (13) значение интеграла будет численно равно значению площади криволинейной трапеции, ограниченной сверху дугой параболы, проходящей через точки: [ х 0 , f 0 )], [ х 1 , f 1 )], [ х 2 , f (х 2 )]

На рисунке сплошной линией изображен график функции f (x ) пунктирной - график многочлена Р 2 (х).

Для получения более точного результата достаточно разбить отрезок интегрирования [а; b ] на четное число (2n ) частей и применить формулу (13) для каждой пары смежных отрезков разбиения:

(14)

Суммируя равенства (14), получим обобщенную формулу Симпсона (парабол):

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

Итак,
.

Выше показали, что
.

Абсолютная погрешность найденного значения не превосходит 0,000005.

Сравнение приближенных значений интеграла , вычисленных по разным формулам, показывает, что наиболее точное значение было получено по обобщенной формуле Симпсона и наименее точное - по формуле прямоугольников.

Погрешность r обобщенной формулы Симпсона можно вычислить по формуле

(16)

где а < ξ< b.

Для абсолютной погрешности обобщенной формулы Симпсона можно получить следующую оценку:

где
(17)

Сравнение точности квадратурных формул.

Выше были приведены оценки абсолютной погрешности квадратурных формул:

для формул прямоугольников: |r|
;

для обобщенной формулы трапеции: |r|
;

для обобщенной формулы Симпсона: |r|
,

где М i =
|f (i) (x)|.

Сопоставление этих оценок позволяет сделать следующие выводы:

    Т.к. производная порядка n+1 от многочлена степениnравна нулю, то получаем точно значение интеграла: по формулетрапеций , если подынтегральная функция линейна,

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

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

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

Для иллюстрации сделанных выводов обратимся к сравнению результатов вычисления интеграла

по различным квадратурным формулам. Для оценки погрешностей вычислим производные функции
.

На отрезке все производные являются монотонными функциями. Абсолютная величина каждой из них достигает своего наибольшего значения при x=0, поэтому М 1 =1, М 2 =2, М 4 =24.

Это позволяет получить при вычислении соответствующие оценки погрешностей:

по формуле прямоугольников r≤0,05;

по формуле трапеций r≤ 0,0017;

по формуле Симпсона r≤ 0,000033.

Сравним полученные результаты, полученные по разным квадратурным формулам со значением ln20,6931472:

по формуле прямоугольников 0,71877;

по формуле трапеций 0,69377;

по формуле Симпсона 0,69315

Видно, что оценки погрешности, как и следовало, ожидать, оказались несколько завышенными.

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

Практические приемы оценки погрешности вычислений по квадратурным формулам.

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

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

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

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

1. Вычисление интеграла
выполняется два раза с шагамиh и 2h . Полученные значения интеграла обозначаются соответственно S h и S 2 h .

2. Если предположить, что на рассматриваемом отрезке [а; b] вторая производная f "(x ) изменяется медленно, то при вычислении интеграла по формуле трапеций можно воспользоваться следующим приближенным выражением для погрешности:

(18)

3. В качестве исправленного (приближенного) значения интеграла можно взять следующее значение:

(19)

    Если предположить, что на рассматриваемом отрезке [а; b] четвертая производная f (4) (х) изменяется медленно, то при вычислении интеграла по формуле Симпсона можно считать, что погрешность приближенно равна

(20)

В качестве исправленного (приближенного) значения интеграла в этом случае можно взять:

(21)

В вычислительной практике часто пользуются также следующим правилом подсчета верных знаков в полученном результате: считают практически верными все совпадающие цифры значений S h иS 2 h .

Приближенное вычисление площадей плоских фигур

Пусть плоская фигура Р ограничена замкнутым контуром С. Выберем систему координат таким образом, чтобы рассматриваемая фигура лежала в пером квадранте. Будем предполагать, что любая прямая, параллельная осиОу, пересекает контур С не более, чем в двух точках. Спроецируем фигуру Р на осьОх ; в проекции получится отрезок[ a ; b ] .

Пусть А – точка фигуры с абсциссой х = а , В – точка фигуры с абсциссойх = b . Точки А и В разбивают контур С на две кривые верхнюю и нижнюю с уравнениями соответственноy = f (x ) иy = g (x ), гдеf (x ), g (x ) – непрерывные на отрезке[ a ; b ] функции. Обозначим черезР площадь фигуры Р. ПлощадьР будет равна разности площадей двух криволинейных трапеций:

аАтВ b иaAhBb ,

т.е. численно равна разности двух интегралов:

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

Разобьем отрезок [а; b ] наn равных частей

0 , х 1 ] , [х 1 , х 2 ], …,[ х п-1 ; х п ]

(а=х 0 , х 1 , …, х п = b ).

Значения подынтегральной функции y = f (x ) - g (x ) будут вычисляться в узлах квадратурной формулы по соотношениям:

y i = f(x i ) - g(x i ) (i = 0, 1, …, п ) .

Очевидно, что

y 0 = f (x 0 ) - g (x 0 ) = 0 и y n = f (x n ) - g (x n ) = 0

Значения y i – длины отрезков ординат в узловых точках, заключенных внутри фигуры Р. Если аналитические выражения функцийf (x ) иg (x ) неизвестны, тоy i можно измерить, пользуясь чертежом.

Общие формулы Ньютона-Котеса

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

I =
,

если на отрезке [а; b ] функция задана таблицей спостоянным шагомh :

x i

x 0

x 1

x 2

x n

y i

y 0

y 1

y 2

y n

Подынтегральную функцию заменим первым интерполяционным многочленом Ньютона и получим:

f (x ) = P n (x ) + R n (x ) (22)

где R n (x ) – остаточный член интерполирования. Интегрируя равенство (22), получим:

отбрасывая второе слагаемое в правой части, получим приближенное равенство

, (23)

погрешность которого определяется формулой:

. (24)

Равенство (23) называют квадратурными формулами Ньютона-Котеса. Из формулы (23) прип=1 получается формула трапеций, а прип =2 – формула Симпсона.

Вычисление интегралов простейшим методом Монте-Карло

Каким образом с помощью кучи камней измерить площадь пруда? Предположим, что пруд расположен в центре поля известной площади А. Бросайте камни в пруд произвольным образом так, чтобы они падали в случайных точках в пределах поля, и считайте количество всплесков при попадании камней в пруд. Эта простая процедура является примером метода Монте-Карло.

Выясним подробнее суть этого метода. Пусть дан прямоугольник высотойН и длинойb - a такой, что функцияf (x ) лежит внутри него. Генерируемп пар случайных чиселx i иy i , удовлетворяющих условиямa <= x i <= b и0 <= y i <= H . Доля точек(x i , y i ) , которые удовлетворяют условиюy i <=f (x i ) , представляет собой оценку отношения интеграла от функцииf (x ) к площади прямоугольника. Отсюда оценкаF n в методе "проб и ошибок" определяется выражением

, (4)

где n s число "всплесков" или точек, лежащих под кривой,п – общее количество точек, а А – площадь прямоугольника.

Другая разновидность метода Монте-Карло основывается на теореме математического анализа, согласно которой определенный интеграл

определяется средним значением подынтегральной функции f (x ) на отрезке[ a ; b ]. Для вычисления этого среднего возьмемx i не с постоянным шагом, а случайным образом и произведемвыборку значенийf (x ) . ОценкаF n одномерного интеграла

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

Численное интегрирование (историческое название: (численная) квадратура ) - вычисление значения определённого интеграла (как правило, приближённое). Под численным интегрированием понимают набор численных методов отыскания значения определённого интеграла.

Численное интегрирование применяется, когда:

В этих двух случаях невозможно вычисление интеграла по формуле Ньютона-Лейбница . Также возможна ситуация, когда вид первообразной настолько сложен, что быстрее вычислить значение интеграла численным методом.

Одномерный случай

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

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

Частным случаем является метод построения интегральных квадратурных формул для равномерных сеток, известный как формулы Котеса . Метод назван в честь Роджера Котса . Основной идеей метода является замена подынтегральной функции каким-либо интерполяционным многочленом . После взятия интеграла можно написать

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

Частными случаями формул Котеса являются: формулы прямоугольников (n=0), формулы трапеций (n=1), формула Симпсона (n=2), формула Ньютона (n=3) и т. д.

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

Пусть требуется определить значение интеграла функции на отрезке . Этот отрезок делится точками на равных отрезков длиной Обозначим через значение функции в точках Далее составляем суммы Каждая из сумм - интегральная сумма для на и поэтому приближённо выражает интеграл

Если заданная функция - положительная и возрастающая, то эта формула выражает площадь ступенчатой фигуры, составленной из «входящих» прямоугольников, также называемая формулой левых прямоугольников, а формула

выражает площадь ступенчатой фигуры, состоящей из «выходящих» прямоугольников, также называемая формулой правых прямоугольников. Чем меньше длина отрезков, на которые делится отрезок , тем точнее значение, вычисляемое по этой формуле, искомого интеграла.

Очевидно, стоит рассчитывать на бо́льшую точность если брать в качестве опорной точки для нахождения высоты точку посередине промежутка. В результате получаем формулу средних прямоугольников:

Учитывая априорно бо́льшую точность последней формулы при том же объеме и характере вычислений её называют формулой прямоугольников

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

Если функцию на каждом из частичных отрезков аппроксимировать прямой, проходящей через конечные значения, то получим метод трапеций.

Площадь трапеции на каждом отрезке:

Погрешность аппроксимации на каждом отрезке:

где

Полная формула трапеций в случае деления всего промежутка интегрирования на отрезки одинаковой длины :

где

Погрешность формулы трапеций:

где

Метод парабол (метод Симпсона)

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

.

Если разбить интервал интегрирования на равных частей, то имеем

Увеличение точности

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

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

При стремлении количества разбиений к бесконечности, оценка интеграла стремится к его истинному значению для аналитических функций для любого численного метода.

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

Метод Гаусса

Описанные выше методы используют фиксированные точки отрезка (концы и середину) и имеют низкий порядок точности (1 - методы правых и левых прямоугольников, 2 - методы средних прямоугольников и трапеций, 3 - метод парабол (Симпсона)). Если мы можем выбирать точки, в которых мы вычисляем значения функции , то можно при том же количестве вычислений подынтегральной функции получить методы более высокого порядка точности. Так для двух (как в методе трапеций) вычислений значений подынтегральной функции, можно получить метод уже не 2-го, а 3-го порядка точности:

.

В общем случае, используя точек, можно получить метод с порядком точности . Значения узлов метода Гаусса по точкам являются корнями полинома Лежандра степени .

Значения узлов метода Гаусса и их весов приводятся в справочниках специальных функций. Наиболее известен метод Гаусса по пяти точкам.

Метод Гаусса-Кронрода

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

,

где - узлы метода Гаусса по точкам, а параметров , , подобраны таким образом, чтобы порядок точности метода был равен .

Тогда для оценки погрешности можно использовать эмпирическую формулу :

,

где - приближённое значение интеграла, полученное методом Гаусса по точкам. Библиотеки gsl и SLATEC для вычисления определённых интегралов содержат подпрограммы, использующие метод Гаусса-Кронрода по 15, 21, 31, 41, 51 и 61 точкам. Библиотека использует метод Гаусса-Кронрода по 15 точкам.

Метод Чебышева

Интегрирование при бесконечных пределах

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

Методы Монте-Карло

Рисунок 3. Численное интегрирование функции методом Монте-Карло

Для определения площади под графиком функции можно использовать следующий стохастический алгоритм:

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

Методы Рунге-Кутты

Метод сплайнов

Многомерный случай

В небольших размерностях можно так же применять квадратурные формулы, основанные на интерполяционных многочленах . Однако в больших размерностях эти методы становятся неприемлемыми из-за быстрого возрастания числа точек сетки и/или сложной границы области. В этом случае применяется метод Монте-Карло . Генерируются случайные точки в нашей области и усредняются значения функции в них. Так же можно использовать смешанный подход - разбить область на несколько частей, в каждой из которых (или только в тех, где интеграл посчитать не удаётся из-за сложной границы) применить метод Монте-Карло .

Литература

  1. Каханер Д., Моулер К., Нэш С. Численные методы и программное обеспечение (пер. с англ.). М.: Мир, 2001, 575 c.

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

Введение

2. Квадратурные формулы

3. Автоматический выбор шага интегрирования

Заключение

Библиографический список


Введение

Цель реферата состоит в изучение и сравнительный анализ методов численного интегрирования функций; реализация этих методов в виде машинных программ на языке высокого уровня и практическое решение задач численного интегрирования на ЭВМ.

При решении инженерных задач часто возникает необходимость в вычислениях значений определенного интеграла вида

Если функция непрерывна на отрезке [a , b ] и ее первообразная может быть определена через известную функцию, то вычисление такого интеграла производится по формуле Ньютона – Лейбница:

.

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

Идея таких методов заключается в следующем. Вместо того, чтобы вычислять интеграл по формуле (1), сначала вычисляют значения функции f (x i ) = y i в некоторых узлах x i Î[a , b ]. Затем выбирается интерполяционный многочлен P (x ), проходящий через полученные точки (x i , y i ), который используется при вычислении приближенного значения интеграла (1):

.

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

, (2)

где - узлы интерполирования, A i – некоторые коэффициенты, R – остаточный член, характеризующий погрешность формулы. Заметим, что формулы вида (2) называют квадратурными формулами.

Геометрический смысл численного интегрирования состоит в вычислении площади криволинейной трапеции, ограниченной графиком функции f (х ),осью абсцисс и двумя прямыми х = а и х = b. Приближенное вычисление площади приводит к отбрасыванию в квадратурных формулах остаточного члена R , характеризующего погрешность метода, на которую дополнительно накладывается вычислительная погрешность.


Методы численного интегрирования

В прикладных исследованиях часто возникает необходимость вычисления значения определённого интеграла

Как известно из курса математики, аналитически вычисление интеграла можно провести не во всех случаях. И даже в том случае, когда удаётся найти аналитический вид этого интеграла, процедура вычисления даёт приближённый результат, поэтому возникает задача приближенного значения этого интеграла.

Суть приближенного вычисления заключается в двух операциях: 1. в выборе конечного числа вместо n; 2. в выборе точки в соответствующем отрезке.

В зависимости от выбора мы получаем различные формулы для вычисления интеграла: Формулы левых и правых прямоугольников (5), (6)

(5)

(6)

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


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

b, a - концы рассматриваемого отрезка.

Для сравнения результатов вычисления вышеизложенными формулами численного интегрирования вычислим 3-мя способами следующий интеграл, разделив отрезок на 6 равных отрезков:

По формуле левых прямоугольников:

По формуле трапеции:

По формуле Симпсона:


А результат полученный аналитически равен

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

Квадратурные формулы

Формулы прямоугольников являются наиболее простыми квадратурными формулами. Разобьем отрезок интегрирования [a, b ] на п равных частей длиной . Заметим, что величину h называют шагом интегрирования. В точках разбиения х 0 = а , х 1 = a + h , ..., x n = b отметим ординаты y 0 , y 1 ,…, y n кривой f (x ), т.е. вычислим у i = f (x i ), x i = a+ ih = x i -1 + h (i = ). На каждом отрезке длиной h построим прямоугольник со сторонами h и y i , где i = , т.е. по значениям ординат, вычисленных в левых концах отрезков. Тогда площадь криволинейной трапеции, определяющую величину интеграла (1), приближенно можно представить в виде суммы площадей прямоугольников (рис. 1). Отсюда получим формулу прямоугольников:


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

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

. (5)

Формулы (3), (4) и (4) называют формулами левых, правых и центральных прямоугольников соответственно.





Рис. 2

Формула трапеций. Здесь на каждом элементарном интервале [x i -1 , x i ] длины h точки с координатами (x i -1 , y i -1) и (x i , y i ) соединяются отрезком (рис. 3). Тогда площадь трапеции, построенной на этом интервале, определяется произведением 0,5h (y i -1 + y i ). Суммируя площади элементарных трапеций для i = получим приближенное значение интеграла.

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

Введение

1. Методы численного интегрирования

2. Квадратурные формулы

3. Автоматический выбор шага интегрирования

Заключение

Библиографический список

Введение

Цель реферата состоит в изучение и сравнительный анализ методов численного интегрирования функций; реализация этих методов в виде машинных программ на языке высокого уровня и практическое решение задач численного интегрирования на ЭВМ.

При решении инженерных задач часто возникает необходимость в вычислениях значений определенного интеграла вида

. (1)

Если функция непрерывна на отрезке [a , b ] и ее первообразная может быть определена через известную функцию, то вычисление такого интеграла производится по формуле Ньютона – Лейбница:

.

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

Идея таких методов заключается в следующем. Вместо того, чтобы вычислять интеграл по формуле (1), сначала вычисляют значения функции f (x i ) = y i в некоторых узлах x i Î[a , b ]. Затем выбирается интерполяционный многочлен P (x ), проходящий через полученные точки (x i , y i ), который используется при вычислении приближенного значения интеграла (1):

.

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

, (2) - узлы интерполирования, A i – некоторые коэффициенты, R – остаточный член, характеризующий погрешность формулы. Заметим, что формулы вида (2) называют квадратурными формулами.

Геометрический смысл численного интегрирования состоит в вычислении площади криволинейной трапеции, ограниченной графиком функции f (х ),осью абсцисс и двумя прямыми х = а и х = b. Приближенное вычисление площади приводит к отбрасыванию в квадратурных формулах остаточного члена R , характеризующего погрешность метода, на которую дополнительно накладывается вычислительная погрешность.

1. Методы численного интегрирования

В прикладных исследованиях часто возникает необходимость вычисления значения определённого интеграла

Как известно из курса математики, аналитически вычисление интеграла можно провести не во всех случаях. И даже в том случае, когда удаётся найти аналитический вид этого интеграла, процедура вычисления даёт приближённый результат, поэтому возникает задача приближенного значения этого интеграла.

Суть приближенного вычисления заключается в двух операциях: 1. в выборе конечного числа вместо n; 2. в выборе точки

в соответствующем отрезке.

В зависимости от выбора

мы получаем различные формулы для вычисления интеграла: Формулы левых и правых прямоугольников (5), (6) (5) (6)

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


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

b, a - концы рассматриваемого отрезка.

Для сравнения результатов вычисления вышеизложенными формулами численного интегрирования вычислим 3-мя способами следующий интеграл, разделив отрезок на 6 равных отрезков: h=

По формуле левых прямоугольников:

По формуле трапеции:

По формуле Симпсона:


А результат полученный аналитически равен

=1

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

2. Квадратурные формулы

Формулы прямоугольников являются наиболее простыми квадратурными формулами. Разобьем отрезок интегрирования [a, b ] на п равных частей длиной

. Заметим, что величину h называют шагом интегрирования. В точках разбиения х 0 = а , х 1 = a + h , ..., x n = b отметим ординаты y 0 , y 1 ,…, y n кривой f (x ), т.е. вычислим у i = f (x i ), x i = a+ ih = x i -1 + h (i = ). На каждом отрезке длиной h построим прямоугольник со сторонами h и y i , где i = , т.е. по значениям ординат, вычисленных в левых концах отрезков. Тогда площадь криволинейной трапеции, определяющую величину интеграла (1), приближенно можно представить в виде суммы площадей прямоугольников (рис. 1). Отсюда получим формулу прямоугольников:
. (3)

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

. (4)

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

. (5)

Формулы (3), (4) и (4) называют формулами левых, правых и центральных прямоугольников соответственно.




Формула Симпсона. Разобьем интервал интегрирования на 2n равных частей длиной

. На каждом отрезке [x i , x i+2 ] подынтегральную функцию f (х ) заменим параболой, проходящей через точки (x i , y i ), (x i +1 , y i +1), (x i +2 , y i +2). Тогда приближенное значение интеграла определяется формулой Симпсона: . (7)

При вычислениях на ЭВМ более удобна следующая формула:


Метод Симпсона - один из наиболее широко известных и применяемых методов численного интегрирования, он дает точные значения интеграла при интегрировании многочленов до третьего порядка включительно.

Формула Ньютона. Приближенное значение интеграла по формуле Ньютона вычисляется следующим образом:

где число участков разбиения кратно трем, т.е. составляет 3n . При разработке программ для ЭВМ удобнее использовать эквивалентную формулу:


Метод Ньютона дает точные значения интеграла при интегрировании многочленов до четвертого порядка включительно.

3. Автоматический выбор шага интегрирования

В результате расчета по формулам (3) - (8) получают приближенное значение интеграла, которое может отличаться от точного на некоторую величину, называемую погрешностью интегрирования. Ошибка определяется формулой остаточного члена R , различной для каждого из методов интегрирования. Если требуется вычислить значение интеграла с погрешностью, не превышающей e, то необходимо выбрать такой шаг интегрирования h , чтобы выполнялось неравенство R (h ) £e. На практике используют автоматический выбор значения h , обеспечивающего достижение заданной погрешности. Сначала вычисляют значение интеграла I (n ), разбивая интервал интегрирования на п участков, затем число участков удваивают и вычисляют интеграл I (2n ). Процесс вычислений продолжают до тех пор, пока не станет справедливым условие.

Страница 1

Кафедра «Высшей математики»
Реферат:

Выполнил: Матвеев Ф.И.
Проверила: Бурлова Л.В.

Улан-Удэ.2002

1.Численные методы интегрирования

2.Вывод формулы Симпсона

3.Геометрическая иллюстрация

4.Выбор шага интегрирования

5.Примеры

1. Численные методы интегрирования
Задача численного интегрирования заключается в вычислении интеграла

посредством ряда значений подынтегральной функции
.

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

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

Численные методы условно можно сгруппировать по способу аппроксимации подынтегральной функции.

Методы Ньютона-Котеса основаны на аппроксимации функции
полиномом степени . Алгоритм этого класса отличается только степенью полинома. Как правило, узлы аппроксимирующего полинома – равноотносящие.

Методы сплайн-интегрирования базируются на аппроксимации функции
сплайном-кусочным полиномом.

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

Методы Монте-Карло используются чаще всего при вычислении кратных интегралов, узлы выбираются случайным образом, ответ носит вероятностный характер.



суммарная погрешность

погрешность усечения

погрешность округления

Независимо от выбранного метода в процессе численного интегрирования необходимо вычислить приближенное значение интеграла и оценить погрешность. Погрешность уменьшается при увеличении n-количества

разбиений отрезка
. Однако при этом возрастает погрешность округления

за счет суммирования значений интегралов, вычисленных на частичных отрезках.

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

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

Проинтегрируем
:

Формула:


и называется формулой Симпсона.

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

Оценим теперь погрешность интегрирования по формуле Симпсона. Будем считать, что у на отрезке
существуют непрерывные производные
. Составим разность

К каждому из этих двух интегралов уже можно применить теорему о среднем, поскольку
непрерывна на
и функция неотрицательна на первом интервале интегрирования и неположительна на втором (то есть не меняет знака на каждом из этих интервалов). Поэтому:

(мы воспользовались теоремой о среднем, поскольку
- непрерывная функция;
).

Дифференцируя
дважды и применяя затем теорему о среднем, получим для
другое выражение:

, где

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

,
.

Если отрезок
интегрирования слишком велик, то его разбивают на
равных частей (полагая
), после чего к каждой паре соседних отрезков
,
,...,
применяют формулу Симпсона, именно:

Запишем формулу Симпсона в общем виде:

(1)

(2)

Погрешность формулы Симпсона - метода четвертого порядка:

,
(3)

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

Например, для функции форма трапеции при
для
дает точный результат
, тогда как по формуле Симпсона получаем

3. Геометрическая иллюстрация


На отрезке
длиной 2h строится парабола, проходящая через три точки
,
. Площадь под параболой, заключенная между осью OX и прямыми
, принимают равной интегралу
.

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

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

(4)

Это формула Симпсона «трех восьмых».

Для произвольного отрезка интегрирования
формула (4) может быть «продолжена»; при этом число частичных отрезков должно быть кратно трем (
точек).

, m=2,3,... (5)

- целая часть

Можно получить формулы Ньютона-Котеса старших порядков:

(6)

- количество отрезков разбиения;

- степень используемого полинома;

- производная -го порядка в точке
;

- шаг разбиения.

В таблице 1 выписаны коэффициенты
. Каждая строка соответствует одному набору промежутков
узлами для построения многочлена k-ой степени. Чтобы воспользоваться этой схемой для большего количества наборов (например, при k=2 и n=6), нужно «продолжить» коэффициенты, а затем сложить их.


Таблица 1:

k

C0

A0

a1

a2

a3

a4

a5

a6

2



1

4

1

1

4

1

1

4

1

1

4

2

2

4

1


Алгоритм оценки погрешности формул трапеции и Симпсона можно записать в виде:
(7),

где - коэффициент, зависящий от метода интегрирования и свойств подынтегральной функции;

h - шаг интегрирования;

p - порядок метода.

Правило Рунге применяют для вычисления погрешности путем двойного просчета интеграла с шагами h и kh.

(8)

(8) - апостериорная оценка. Тогда Iуточн.= +Ro (9),
уточненное значение интеграла
.

Если порядок метода неизвестен, необходимо вычислить I в третий раз с шагом
, то есть:

из системы трех уравнений:

с неизвестными I,А и p получаем:

(10)

Из (10) следует
(11)

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

Аналогично, для интегралов, вычисленных по формуле с шагами и
, справедливы соотношения:

,

(13)

4. Выбор шага интегрирования
Для выбора шага интегрирования можно воспользоваться выражением остаточного члена. Возьмем, например, остаточный член формулы Симпсона:

Если 

, то 

.

По заданной точности  метода интегрирования из последнего неравенства определяем подходящий шаг.

,
.

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

Разберем один из таких приемов. Пусть

,

где - приближенное значение интеграла с шагом . Уменьшим шаг в два раза, разбив отрезок
на две равные части
и
(
).

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

Отсюда можно сделать такой вывод: если
, то есть если
,
, а - требуемая точность, то шаг подходит для вычисления интеграла с достаточной точностью. Если же
, то расчет повторяют с шагом и затем сравнивают
и
и т.д. Это правило называется правилом Рунге.

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

При выводе правила Рунге вы существенно пользовались предположением, что
. Если имеется только таблица значений , то проверку
«на постоянство» можно сделать непосредственно по таблице Дальнейшее развитие приведенных алгоритмов позволяет перейти к адаптивным алгоритмам, в которых за счет выбора различного шага интегрирования в разных частях отрезка интегрирования в зависимости от свойств
уменьшается количество вычислений подынтегральной функции.

Другая схема уточнения значений интеграла - процесс Эйтнена. Производится вычисление интеграла с шагами
, причем
. Вычисление значений . Тогда
(14).

За меру точности метода Симпсона принимают величину:

5. Примеры
Пример 1. Вычислить интеграл
по формуле Симпсона, если
задана таблицей. Оценить погрешность.

Таблица 3.




0

0.1

0.2

0.3

0.4

0.5

0.6

0.7

0.8




1

0.995

0.98

0.955

0.921

0.878

0.825

0.765

0.697

Решение: Вычислим по формуле (1) при
и
интеграл .

По правилу Рунге получаем
Принимаем .


Пример 2. Вычислить интеграл
.

Решение: Имеем
. Отсюда h=
=0.1. Результаты вычислений приведены в таблице 4.


Таблица 4.

Вычисление интеграла по формуле Симпсона


i









0

0

y0=1,00000

1

0.1

0,90909

2

0.2

0,83333

3

0.3

0,76923

4

0.4

0,71429

5

0.5

0,66667

6

0.6

0,62500

7

0.7

0,58824

8

0.8

0,55556

9

0,9

0,52632

10

1,0

0,50000=yn



3,45955(1)

2,72818(2)

По формуле Симпсона получим:

Подсчитаем погрешность полученного результата. Полная погрешность складывается из погрешностей действий и остаточного члена . Очевидно:-0,289687


4

2,35

-0,70271

-0,299026

2,4

-0,73739

-0,307246

2

2,45

-0,77023

-0,314380

2,5

-0,80114

-0,320465

4

2,55

-0,83005

-0,325510

2,6

-0,85689

-0,329573

2

2,65

-0,88158

-0,332672

2,7

-0,90407

-0,334841

4

2,75

-0,92430

-0,336109

 3.